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Improved  Ultrasonic  Imaging  of  the  Breast 
William  F.  Walker,  Ph.D. 

Technical  Abstract: 

Ultrasonic  imaging  has  become  an  increasingly  important  tool  for  the  diagnosis  of  breast  cancer.  While 
x  ray  mammography  remains  the  standard  for  screening,  ultrasound  is  widely  used  to  differentiate  fluid  filled 
cysts  from  solid  masses,  to  guide  invasive  biopsy  procedures,  and  even  to  differentiate  between  benign  and 
malignant  lesions.  Ultrasound  is  an  attractive  choice  for  these  applications  because  it  exposes  the  patient  to  no 
ionizing  radiation,  requires  no  uncomfortable  breast  compression,  produces  images  in  real-time,  and  can 
successfully  image  young  women  with  radiographically  dense  breasts.  Although  it  offers  many  advantages, 
ultrasound  is  limited  by  the  fact  that  certain  soft  tissues  may  appear  only  subtly  different  and  that 
microcalcifications  (MCs)  are  not  reliably  visualized.  The  inability  to  reliably  image  MCs  is  particularly 
troublesome  because  these  small  calcium  crystals  are  an  important  mammographic  indicator  of  breast  cancer. 

This  grant  would  support  career  development  by  reducing  teaching  load  to  provide  time  for  expanded 
research  efforts  and  clinical  training.  If  funded,  my  teaching  load  would  fall  from  three  courses  per  academic 
year  to  one  course  per  academic  year.  I  have  been  supported  by  the  University  of  Virginia  by  provision  of 
excellent  laboratory  space  within  the  Biomedical  Engineering  Department,  Medical  Research  Building  4,  and 
occasionally  in  the  Department  of  Radiology.  This  laboratory  space  is  well  equipped  with  test  equipment,  an 
ophthalmic  ultrasound  system,  a  custom  ultrasound  system,  and  a  GE  Logiq  700  ultrasonic  imaging  system  with 
research  interface.  The  custom  system  has  precise  timing  control,  making  it  an  excellent  platform  for 
development  of  new  imaging  and  flow  estimation  algorithms.  The  GE  system  is  coupled  to  an  extended 
research  interface  which  allows  acquisition  of  up  to  32  MB  of  raw  ultrasound  echo  data.  This  system  is  also 
supported  by  custom  research  software,  allowing  control  of  aperture  geometry  and  a  variety  of  other  system 
parameters. 

The  main  goals  of  this  grant  are: 

(1)  Research  in  angular  scatter  imaging  using  the  translating  apertures  algorithm. 

(2)  Clinical  training  in  breast  cancer  imaging  and  treatment. 

The  first  aim  will  apply  the  Translating  Apertures  Algorithm  (TAA)  to  form  images  of  the  angular 
scatter  parameters  of  tissue.  Since  the  TAA  observes  the  same  speckle  pattern  with  angle  of  interrogation  (for 
omnidirectional  scatterers),  it  eliminates  the  statistical  variability  which  plagued  earlier  techniques.  This 
property  makes  it  possible  for  the  TAA  to  acquire  independent,  statistically  reliable  angular  scatter  profiles  at 
every  location  in  tissue.  The  specific  aims  of  this  research  are: 

1.  Determine  bias  and  variance  of  angular  scatter  measurements  performed  with  the  TAA. 

2.  Form  novel  images  using  angular  scatter  data  from  multiple  angles. 

3.  Implement  angular  scatter  imaging  in  combination  with  spatial  compounding. 

Each  of  these  aims  will  be  performed  using  experimental  phantom  data,  experimental  human  tissue  data,  and 
computer  simulations.  All  experiments  will  be  performed  on  a  GE  Logiq  700  imaging  system  with  a  special 
research  interface.  The  images  obtained  in  aims  2  and  3  will  be  formed  offline  using  experimental  data  obtained 
from  the  GE  Logic  700  system  described  above. 

The  techniques  developed  under  this  grant  will  likely  improve  the  visibility  of  MCs  and  should  also 
enhance  the  contrast  of  soft  tissue  lesions.  Since  tumors  are  know  to  have  different  concentrations  of 
extracellular  matrix  proteins  they  would  be  expected  to  exhibit  angular  scatter  responses  different  from  that  of 
normal  tissue.  Since  MCs  have  a  much  higher  mass  density  that  their  surrounding  tissue,  they  would  be 
expected  to  exhibit  a  much  greater  angular  scatter  variation  than  soft  tissues. 


We  have  already  begun  development  of  one  method  of  angular  scatter  imaging.  In  this  method  we  form 
two  separate  images  of  tissue;  one  which  highlights  the  component  of  scattering  which  is  uniform  with  angle  (c- 
weighted  image),  and  a  second  which  highlights  the  scattering  component  which  varies  with  angle  (d-weighted 
image).  Simulation  results  indicate  that  c-  and  d-weighted  images  may  offer  significant  new  information  about 
soft  tissue,  and  that  they  will  almost  certainly  improve  the  detectibility  of  MCs  with  ultrasound. 

The  techniques  developed  here  should  allow  detection  of  previously  invisible  tumors,  especially  in 
women  with  radiodense  breasts.  Furthermore  they  should  improve  the  accuracy  of  image  guided  biopsy 
procedures  by  enhancing  the  visibility  of  MCs. 

I  believe  that  it  is  my  job  to  develop  new  technologies  which  can  be  implemented  to  improve  patient 
care.  While  my  knowledge  of  engineering  continues  to  evolve,  I  have  become  increasingly  aware  of  limitations 
in  my  knowledge  of  clinical  medicine.  If  this  proposal  is  funded  I  will  invest  a  significant  amount  of  time  to 
educate  myself  in  the  clinical  methods  of  breast  cancer  diagnosis  and  treatment.  I  will  observe  breast  imaging, 
biopsy,  and  surgical  procedures  to  gain  first  hand  knowledge  of  the  strengths  and  weaknesses  of  existing 
methods.  I  believe  that  these  observations  will  illuminate  new  directions  of  research. 


Narrative: 

The  main  goal  of  this  project  is  the  validation  of  angular  scatter  ultrasonic  imaging  as  a  useful  method  for 
improving  the  detection  of  microcalcifications  and  the  differentiation  of  benign  and  malignant  breast  tissues.  In 
the  first  year  of  support  we  successfully  implemented  the  Translating  Apertures  Algorithm  (TAA)  on  the  GE 
Logiq  700MR  imaging  system  located  in  the  Pi’s  laboratory.  We  also  showed  that  angular  scatter  imaging 
obtains  information  about  both  phantoms  and  human  tissues  that  is  not  available  through  existing  methods.  In 
phantom  experiments  we  showed  that  angular  scatter  imaging  improved  the  contrast  of  microcalcification 
mimicking  phantom  materials  by  as  much  as  a  factor  of  5. 

Unfortunately  our  early  in  vivo  work  did  not  show  the  impressive  results  we  expected.  In  the  past  year  we  have 
identified  a  number  of  problems  that  limited  in  vivo  performance.  Our  spatial  resolution  has  been  limited 
because  a  relative  small  available  aperture  requires  us  to  employ  a  small  imaging  aperture,  and  thus  achieve 
poor  spatial  resolution.  We  attempted  to  circumvent  this  problem  by  increasing  the  effective  aperture  size,  but 
found  that  as  the  aperture  size  increased  (increasing  spatial  resolution)  our  ability  to  resolve  variations  in 
angular  scatter  was  lost.  Thus  we  identified  a  fundamental  limit  of  angular  scatter  imaging  with  conventional 
methods:  It  is  impossible  to  simultaneously  achieve  high  spatial  resolution  and  high  angular  resolution.  In 
addition,  our  imaging  system  was  not  able  to  apply  apodization  so  large  sidelobes  greatly  degraded  image 
contrast.  Finally,  throughout  most  of  this  second  year  we  were  plagued  by  phase  errors  between  different 
scattering  angles  that  were  induced  by  hidden  parameters  within  the  GE  system. 

Much  of  the  second  year  of  this  grant  has  been  devoted  to  correcting  for  GE  system  effects  and  to  developing 
techniques  to  simultaneously  improve  spatial  resolution  and  angular  resolution.  The  GE  system  effects  were 
finally  completely  accounted  for  in  August  after  extensive  collaboration  with  K.  Wayne  Rigby  at  the  GE  Global 
Research  Center.  Dr.  Rigby  worked  with  the  PI  to  identify  and  correct  for  a  number  of  hidden  system  effects.  In 
addition  to  the  effects  that  Dr.  Rigby  identified,  there  were  a  few  other  effects  that  could  only  be  corrected  by 
empirical  study.  We  have  now  completed  a  series  of  control  experiments  and  confirmed  that  the  there  are  no 
significant  residual  system  effects. 

To  simultaneously  optimize  spatial  resolution,  image  contrast,  and  angular  resolution  we  have  developed  a 
novel  synthetic  aperture  imaging  method  for  angular  scatter  imaging.  We  have  implemented  this  technique  on 
the  GE  system  and  performed  initial  experiments  indicating  that  we  are  in  fact  able  to  obtain  excellent  spatial 
and  angular  resolution.  Unfortunately  this  approach  has  a  requirement  on  spatial  array  sampling  which  is  twice 
as  restrictive  as  that  seen  for  conventional  beamforming.  We  have  followed  two  strategies  to  address  this 
problem.  In  the  first,  we  modify  our  transmit  pulse  to  utilize  lower  frequency  components  and  thereby  reduce 
the  magnitude  of  grating  lobes  that  would  result  from  imaging  at  the  typical  frequencies  for  our  arrays.  In  our 
second  approach  we  slightly  reduce  the  angular  resolution  to  further  reduce  grating  lobe  magnitudes.  So  far 
these  methods  have  combined  to  reduce  grating  lobe  magnitudes  by  roughly  40  dB,  or  a  factor  of  100.  Our  final 
testing  and  implementation  of  this  technique  will  be  complete  in  the  next  few  weeks.  The  developed  synthetic 
aperture  methods  will  allow  us  to  form  control  images  that  rival  those  of  the  best  ultrasound  imaging  systems 
available  today.  However,  we  will  also  be  able  to  obtain  angular  scatter  images  with  similar  contrast  and 
resolution.  Once  this  method  is  functional  we  will  repeat  our  phantom  experiments  and  continue  in  vivo 
experiments. 

In  the  coming  year  we  will  continue  our  experimental  work  and  test  many  of  our  signal  processing 
methods  experimentally.  We  will  modify  the  operation  of  the  GE  system  to  enhance  image  resolution  while 
maintaining  angular  scatter  resolution.  We  will  also  develop  tools  for  angular  scatter  imaging  on  the  Philips 


SONOS  5500  imaging  system  now  in  our  laboratory.  The  development  of  these  tools  is  being  led  by  an 
Electrical  Engineering  graduate  student,  Greg  Yukl,  who  is  supported  by  an  NSF  major  research 
instrumentation  grant  for  the  development  of  such  tools.  We  believe  that  the  greater  flexibility  and  signal 
fidelity  of  the  SONOS  system  will  improve  angular  scatter  results  and  allow  more  general  experimentation. 

During  2003  we  have  initiated  some  new  research  directions.  Jake  Mann,  a  graduate  student  under  the  Pi’s 
direction,  completed  an  M.S.  thesis  titled  “A  Constrained  Adaptive  Beamformer  for  Medical  Ultrasound.”  This 
thesis  described  a  number  of  new  methods  that  have  great  potential  to  improve  the  contrast  and  resolution  of 
breast  ultrasound.  Mr.  Mann  presented  early  results  at  an  international  meeting  and  we  have  since  received 
notice  from  the  Army  BCRP  that  an  Idea  Grant  proposal  based  on  this  work  has  been  recommended  for 
funding.  In  addition  to  the  constrained  adaptive  beamforming  work  we  have  also  initiated  the  research  on 
synthetic  aperture  methods  descried  above. 

The  past  year  has  seen  a  significant  change  in  the  makeup  of  our  research  team.  Throughout  the  lifetime  of  this 
project  laboratory  work  has  been  performed  by  doctoral  candidate  M.  Jason  McAllister.  Unfortunately  Mr. 
McAllister’s  progress  has  not  proceeded  at  the  rate  that  the  PI  expected  or  required.  As  a  result,  at  the  end  of 
2002,  Mr.  McAllister  was  presented  with  a  detailed  timeline  describing  expected  progress.  When  he  deviated 
failed  to  meet  this  timeline  in  Spring  of  2003  his  financial  support  was  removed.  This  was  not  a  step  I  took 
lightly  or  relished,  but  I  felt  it  was  necessary  to  move  the  project  along.  I  have  since  recruited  a  new  student  to 
the  project,  Mr.  Drake  Guenther.  Mr.  Guenther  has  made  excellent  progress  this  fall,  leading  implementation 
and  testing  of  the  synthetic  aperture  methods  described  above. 


Key  Research  Accomplishments: 


a)  Determine  bias  and  variance  of  angular  scatter  measurements  performed  with  the  TAA. 

Simulate  bias  and  variance  Months  1-18 

Underway,  initial  results  described  in  the  following  publications: 

“A  Novel  Beamformer  Design  Method  for  Medical  Ultrasound:  Part  II:  Simulation  Results” 
“Novel  Aperture  Design  Method  for  Improved  Depth  of  Field  in  Ultrasound  Imag.  ” 

“A  Novel  Aperture  Design  Method  in  Ultrasound  Imaging” 

“Minimum  Sum  Squared  Error  (MSSE)  Beamformer  Design  Techn.:  Initial  Results” 
“Angular  Scatter  Imaging:  Clinical  Results  and  Novel  Processing  Methods” 

Measure  bias  and  variance  in  phantoms  Months  7-18 

Underway,  initial  results  described  in  the  following  publications: 

“Eval.  Trans.  Apertures  Based  Angular  Scatter  Imaging  on  a  Clin.  Imaging  System” 
“Angular  Scatter  Imaging:  Clinical  Results  and  Novel  Processing  Methods” 

Develop  and  test  methods  for  angle  dependent  weightings 

to  compensate  for  limited  element  angular  response  Months  1 3-24 

Underway,  initial  results  described  in  the  following  publications: 

“A  Novel  Beamformer  Design  Method  for  Medical  Ultrasound:  Part  I:  Theory” 

“A  Novel  Beamformer  Design  Method  for  Medical  Ultrasound:  Part  II:  Simulation  Results” 
“Novel  Aperture  Design  Method  for  Improved  Depth  of  Field  in  Ultrasound  Imag.” 

“A  Novel  Aperture  Design  Method  in  Ultrasound  Imaging” 

“ Minimum  Sum  Squared  Error  (MSSE)  Beamformer  Design  Techn.:  Initial  Results” 
“Angular  Scatter  Imaging:  Clinical  Results  and  Novel  Processing  Methods” 

Develop  and  test  apodization  schemes  to  compensate  for 

apparent  apodization  due  to  element  angular  response  Months  19-36 
Underway,  initial  results  described  in  the  following  publications: 

“A  Novel  Beamformer  Design  Method  for  Medical  Ultrasound:  Part  I:  Theory” 

“A  Novel  Beamformer  Design  Method  for  Medical  Ultrasound :  Part  II:  Simulation  Results” 
“Novel  Aperture  Design  Method  for  Improved  Depth  of  Field  in  Ultrasound  Imag.  ” 

“A  Novel  Aperture  Design  Method  in  Ultrasound  Imaging  ” 

“Minimum  Sum  Squared  Error  (MSSE)  Beamformer  Design  Techn.:  Initial  Results” 
“Angular  Scatter  Imaging:  Clinical  Results  and  Novel  Processing  Methods” 


Develop  expressions  relating  correlation  to  variance  Months  1-12 

We  are  no  longer  working  towards  this  original  goal.  Rather  than  comparing  correlation 
and  variance  we  have  focused  the  relationship  between  correlation  and  sum  squared  error. 
Sum  squared  error  provides  a  more  useful  metric  as  it  can  form  the  basis  of  algorithms  to 
improve  system  performance  and  contains  the  effects  of  both  bias  and  variance.  In  the 
following  initial  papers  we  have  considered  the  impact  of  improved  sum  squared  error  on 
correlation: 

“A  Novel  Beamformer  Design  Method  for  Medical  Ultrasound:  Part  II:  Simulation  Results” 
“Novel  Aperture  Design  Method  for  Improved  Depth  of  Field  in  Ultrasound  Imag.  ” 

“A  Novel  Aperture  Design  Method  in  Ultrasound  Imaging  ” 

“Minimum  Sum  Squared  Error  (MSSE)  Beamformer  Design  Techn.:  Initial  Results” 

“Angular  Scatter  Imaging:  Clinical  Results  and  Novel  Processing  Methods” 

b)  Form  novel  images  using  angular  scatter  data  from  multiple  angles. 

Simulate  angular  scatter  data 

(using  Rayleigh,  Faran,  and  other  models)  Months  13-48 

Although  we  have  initiated  work  on  this  aim,  we  are  not  yet  ready  to  report  results. 

Acquire  data  from  tissue  mimicking  phantoms  Months  7-30 

Underway,  initial  results  described  in  the  following  publications: 

“Eval.  Trans.  Apertures  Based  Angular  Scatter  Imaging  on  a  Clin.  Imaging  System” 

“A  Constrained  Adaptive  Beamformer  for  Medical  Ultrasound:  Initial  Results” 

“Angular  Scatter  Imaging:  Clinical  Results  and  Novel  Processing  Methods” 

Acquire  data  from  human  breast  tissue  Months  13-36 

We  are  currently  preparing  to  start  work  on  this  aim. 

Form  parametric  images  using  polynomial  fits  Months  16-39 

Although  we  have  initiated  work  on  this  aim,  we  are  not  yet  ready  to  report  results. 

Adapt  Haider’s  method  to  angular  scatter  imaging  Months  19-42 

We  have  not  yet  initiated  work  on  this  aim. 

Image  variance  of  angular  scatter  Months  16-39 

Although  we  have  initiated  work  on  this  aim,  we  are  not  yet  ready  to  report  results. 


Develop  other  imaging  methods  Months  22-48 

Although  we  have  initiated  work  on  this  aim,  we  are  not  yet  ready  to  report  results. 


Test  methods  to  improve  depth  of  field  Months  13-48 

Underway,  initial  results  described  in  the  following  publications: 

“A  Novel  Beamformer  Design  Method  for  Medical  Ultrasound:  Part  II:  Simulation  Results” 
“Eval.  Trans.  Apertures  Based  Angular  Scatter  Imaging  on  a  Clin.  Imaging  System” 

“Novel  Aperture  Design  Method  for  Improved  Depth  of  Field  in  Ultrasound  Imag.  ” 

“A  Novel  Aperture  Design  Method  in  Ultrasound  Imaging” 

“Minimum  Sum  Squared  Error  (MSSE)  Beamformer  Design  Techn.:  Initial  Results” 
“Angular  Scatter  Imaging:  Clinical  Results  and  Novel  Processing  Methods” 

Develop  and  test  methods  to  detect  specular  reflectors  Months  1-24 

We  have  not  yet  initiated  work  on  this  aim. 

c)  Implement  angular  scatter  imaging  in  combination  with  spatial  compounding. 

Implement  TAA  and  spatial  compounding  on  the  GE  Months  25-36 

We  have  not  yet  initiated  work  on  this  aim. 

Develop  methods  to  reduce  the  impact  of  limited  element 

angular  response  Months  3 1  -48 

Underway,  initial  results  described  in  the  following  publications: 

“A  Novel  Beamformer  Design  Method  for  Medical  Ultrasound:  Part  I:  Theory” 

“A  Novel  Beamformer  Design  Method  for  Medical  Ultrasound:  Part  II:  Simulation  Results” 
“Novel  Aperture  Design  Method  for  Improved  Depth  of  Field  in  Ultrasound  Imag.  ” 

“A  Novel  Aperture  Design  Method  in  Ultrasound  Imaging” 

“Minimum  Sum  Squared  Error  (MSSE)  Beamformer  Design  Techn.:  Initial  Results” 

“A  Constrained  Adaptive  Beamformer  for  Medical  Ultrasound:  Initial  Results  ” 

“Angular  Scatter  Imaging:  Clinical  Results  and  Novel  Processing  Methods” 


Test  compounding  with  the  TAA  on  phantoms 

We  have  not  yet  initiated  work  on  this  aim. 

Months  31-48 

Test  compounding  with  the  TAA  on  tissues 

We  have  not  yet  initiated  work  on  this  aim. 

Months  37-48 

Reportable  Outcomes: 

(Note  that  student  authors  are  underlined.) 


Refereed  Publications: 

Ranganathan.  K.  and  W.  F.  Walker,  “A  Novel  Beamformer  Design  Method  for  Medical  Ultrasound:  Part  I: 
Theory,”  IEEE  Trans.  Ultrason.  Ferroelec.  Freq.  Contr.,  IEEE  Trans.  Ultrason.  Ferroelec.  Freq.  Contr,  vol. 
50,  no.  1,  pp.  15-24, 2003. 

Ranganathan.  K.  and  W.  F.  Walker,  “A  Novel  Beamformer  Design  Method  for  Medical  Ultrasound:  Part  II: 
Simulation  Results,”  IEEE  Trans.  Ultrason.  Ferroelec.  Freq.  Contr.,  vol.  50,  no.  1,  pp.  25-39,  2003. 

Conference  Presentations  with  Paper: 

Mann.  J.  A.  and  W.  F.  Walker,  “Constrained  Adaptive  Beamforming:  Point  and  Contrast  Resolution,”  2003 
SPIE  Medical  Imaging  Symposium. 

Walker,  W.  F.  and  M.  J.  McAllister.  “Angular  Scatter  Imaging  in  Medical  Ultrasound,”  presented  at  the 
36th  Asilomar  Conference  on  Signals,  Systems,  and  Computers. 

McAllister.  M.  J.  and  W.  F.  Walker,  “Angular  Scatter  Imaging:  Clinical  Results  and  Novel  Processing 
Methods,”  Proc.  IEEE  Ultrason.  Symp.,  pp.  ***-***,  2002. 

Mann.  J.  A.  and  W.  F.  Walker,  “A  Constrained  Adaptive  Beamformer  for  Medical  Ultrasound:  Initial 
Results,”  Proc.  IEEE  Ultrason.  Symp.,  pp.  ***-***,  2002. 

Ranganathan.  K..  and  W.  F.  Walker,  “The  Minimum  Sum  Squared  Error  (MSSE)  Beamformer  Design 
Technique:  Initial  Results,”  Proc.  IEEE  Ultrason.  Symp.,  pp.  ***-***,  2002. 

Ranganathan.  K.  and  W.  F.  Walker,  “A  Novel  Aperture  Design  Method  in  Ultrasonic  Imaging,”  Proc.  SPIE 
vol.  4687,  pp.  13-24,  Medical  Imaging  2002:  Ultrasonic  Imaging  and  Signal  Processing,  2002. 

McAllister.  M.  J..  K.  W.  Rigby  and  W.  F.  Walker,  “Evaluation  of  Translating  Apertures  Based  Angular 
Scatter  Imaging  on  a  Clinical  Imaging  System,”  Proc.  IEEE  Ultrason.  Symp.,  pp.  1555-8, 2001. 

Ranganathan.  K.  and  W.  F.  Walker,  “A  Novel  Aperture  Design  Method  for  Improved  Depth  of  Field  in 
Ultrasound  Imaging,”  Proc.  IEEE  Ultrason.  Symp.,  pp.  1543-6,  2001. 

Conference  Presentations  with  Abstract: 

Walker,  W.F.,  M.J.  McAllister.  K.W.  Rigby,  “In  Vivo  Imaging  of  Ultrasonic  Angular  Scatter,”  U.S.  Dept, 
of  Defense  Breast  Cancer  Research  Program  Era  of  Hope  2002  Meeting,  September  2002. 

Presentations: 

Oct.  2003  Parent’s  Weekend,  UVA  SEAS,  Charlottesville,  VA“Ultrasound:  The  Basics  and  the 

Future” 

June,  2003  Alumni  Weekend,  UVA  SEAS,  Charlottesville,  VA  “Ultrasound:  The  Basics  and  the 

Future” 

Oct.,  2002  Philips  Research  Center,  Paris  France  “Medical  Ultrasound  Research  at  UVA” 


Mar.,  2002  Phillips  Medical  Sys.,  Bothell,  WA 

“Radiation  Force,  Angular  Scatter,  &  Novel  Beamforming” 

Feb.,  2002  IBM  Watson  Research  Center,  White  Plains,  NY 

“Compute  &  Storage  Intensive  Ultrasound  Research” 

Patents: 

“Constrained  Adaptive  Beamformer  for  Medical  Ultrasound,”  J.  A.  Mann  and  W.  F.  Walker,  Provisional 
patent  Application  filed  October  7, 2002. 

“An  Improved  Beamforming  Method  for  Ultrasound  Imaging,”  K.  Ranganathan  and  W.  F.  Walker, 
Provisional  Patent  filed  October  5, 2001. 

New  Funding  Applications: 

“Angular  Scatter  Ultrasound  Imaging  in  the  Breast”,  the  National  Institutes  of  Health,  $1,794,454,  funding 
requested  from  December  1,  2003  to  November  30,  2008  (20%  effort).  (Proposal  was  scored  at  186  but  not 
funded.  We  will  resubmit  in  2004.) 

“Constrained  Adaptive  Beamforming  for  Improved  Contrast  in  Breast  Ultrasound,”  Congressionally 
Directed  Medical  Research  Program,  $429,488,  funding  requested  from  January  1,  2004  to  December  31, 
2006  (10%  effort).  (Proposal  was  scored  at  1.5  and  recommended  for  funding.) 

“Shear  Elastography  for  Breast  Cancer  Detection  and  Diagnosis,”  Congressionally  Directed  Medical 
Research  Program,  $442,745,  funding  requested  from  January  1,  2004  to  December  31,  2006  (10%  effort). 
(Proposal  was  scored  at  2.6  but  not  recommended  for  funding.) 

“Ultrasonic  Angular  Scatter  Imaging  Using  Synthetic  Aperture  Methods  for  Breast  Imaging,” 
Congressionally  Directed  Medical  Research  Program,  $442,745,  funding  requested  from  January  1,  2004  to 
December  31, 2006  (10%  effort).  (Proposal  was  scored  at  1.5  but  not  recommended  for  funding.) 


Conclusions: 

Angular  scatter  imaging  with  ultrasound  shows  more  promise  now  than  it  did  at  the  initiation  of 
this  grant.  Our  experimental  results  indicate  that  angular  scatter  variations  are  significant  in  both 
phantoms  and  human  tissues,  and  that  they  offer  a  source  of  contrast  that  is  distinct  from  traditional 
backscatter.  Unfortunately,  it  appears  that  the  application  of  conventional  beamforming  methods  may 
obscure  angular  scatter  information  by  effectively  trading  angular  resolution  for  spatial  resolution.  To 
circumvent  this  tradeoff  we  have  invented  novel  synthetic  aperture  methods  that  achieve  both  high 
spatial  resolution  and  high  angular  resolution.  We  have  implemented  these  methods  and  are  in  early 
stages  of  testing  them  experimentally.  In  the  coming  year  we  will  complete  these  tests,  expand  to 
phantom  testing,  and  extend  in  vivo  testing. 
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A  Novel  Beamformer  Design  Method  for 
Medical  Ultrasound.  Part  I:  Theory 

Karthik  Ranganathan  and  William  F.  Walker.  Member,  IEEE 


Abstract — The  design  of  transmit  and  receive  aperture 
weightings  is  a  critical  step  in  the  development  of  ultra¬ 
sound  imaging  systems.  Current  design  methods  are  gener¬ 
ally  iterative,  and  consequently  time  consuming  and  inex¬ 
act.  We  describe  a  new  and  general  ultrasound  beamformer 
design  method,  the  minimum  sum  squared  error  (MSSE) 
technique.  The  MSSE  technique  enables  aperture  design  for 
arbitrary  beam  patterns  (within  fundamental  limitations 
imposed  by  diffraction).  It  uses  a  linear  algebra  formula¬ 
tion  to  describe  the  system  point  spread  function  (psf)  as 
a  function  of  the  aperture  weightings.  The  sum  squared  er¬ 
ror  (SSE)  between  the  system  psf  and  the  desired  or  goal 
psf  is  minimized,  yielding  the  optimal  aperture  weightings. 
We  present  detailed  analysis  for  continuous  wave  (CW)  and 
broadband  systems.  We  also  discuss  several  possible  appli¬ 
cations  of  the  technique,  such  as  the  design  of  aperture 
weightings  that  improve  the  system  depth  of  field,  generate 
limited  diffraction  transmit  beams,  and  improve  the  corre¬ 
lation  depth  of  field  in  translated  aperture  system  geome¬ 
tries.  Simulation  results  are  presented  In  an  accompanying 
paper. 


3.  Introduction 

ULTRASONIC  imaging  is  an  important  medical  diagnos¬ 
tic  tool  that  entails  four  critical  steps.  An  ultrasound 
waveform  is  first  generated  and  transmitted  by  exciting 
a  piezoelectric  transducer  with  a  suitable  electric  signal; 
usually  a  finite  duration  pulse.  This  transmitted  ultra¬ 
sound  beam  propagates  through  tissue,  undergoing  diffrac¬ 
tion  and  attenuation  as  well  as  scattering  and  reflection  at 
tissue  interfaces.  Next,  the  reflected  or  scattered  echoes 
propagate  back  to  the  transducer  where  they  are  received 
and  converted  to  electric  signals.  These  received  echoes 
are  then  processed  and  mapped  to  form  an  image.  Simple 
ultrasound  systems  use  a  single  large  piezoelectric  trans¬ 
ducer  to  generate  and  receive  ultrasound.  However,  state 
of  the  art  systems  use  phased  array  technology  similar  to 
that  used  in  contemporary  radar  and  sonar.  These  sys¬ 
tems  use  transducer  arrays  comprising  a  multitude  of  small 
transducer  elements.  Subsets  of  these  elements  form  the 
active  transducer  aperture  that  is  used  for  transmission 
and  reception.  The  prime  reason  for  the  use  of  arrays  is 
the  fact  that  the  independent  operation  of  each  element 

Manuscript  received  May  16,  2002*,  accepted  August  28,  2002.  This 
work  was  supported  by  Susan  G.  Komen  Breast  Cancer  Foundation 
Imaging  Grant  No.  99-3021  and  United  States  Army  Congrcssionally 
Directed  Medical  Research  Program  Grant  No.  DAMD  17-01-10443. 
Inspiration  from  this  work  stems  from  National  Science  Foundation 
Major  Research  Instrumentation  Grant  0079639. 

The  authors  are  with  the  Department  of  Biomedical  En¬ 
gineering,  University  of  Virginia,  Charlottesville,  VA  (e-mail: 
kr6u@virginia.edu). 


enables  more  control  over  the  transmission  and  reception 
processes.  On  transmit,  the  application  of  distinct  time 
delays  [1]  to  the  pulses  used  to  excite  each  element  fo¬ 
cuses  the  transmitted  ultrasound  beam  to  a  specific  point. 
These  delays  are  calculated  using  elementary  geometry  to 
compensate  for  differences  between  the  path  lengths  from 
each  element  to  the  point  of  interest.  In  addition  to  these 
time  delays,  magnitude  weights  [1]  may  be  applied  to  the 
elements  (apodization)  to  control  the  shape  of  the  ultra¬ 
sonic  beam.  These  weights  are  usually  determined  itera¬ 
tively  or  adapted  from  windowing  functions  described  in 
signal  processing  literature.  These  magnitude  and  phase 
weights  play  a  critical  role  in  determining  system  resolu¬ 
tion  and  contrast.  More  advanced  medical  ultrasound  sys¬ 
tems  enable  not  only  transmit  focusing  and  apodization, 
but  also  the  generation  of  arbitrary  transmit  waveforms 
for  each  element. 

On  the  receive  side,  state  of  the  art  systems  change  their 
focus  dynamically  as  a  function  of  time  after  transmit  [1], 
This  means  that  the  system  is  ideally  focused  at  the  point 
of  origin  of  the  echoes  received  at  any  time.  In  addition  to 
dynamic  focusing,  these  systems  also  implement  dynamic 
apodization  (l]  to  maintain  a  constant  system  f-number 
within  the  physical  constraints  of  the  transducer.  The  f- 
number  is  the  ratio  of  the  range  being  interrogated  to  the 
aperture  size.  Keeping  the  f-number  constant  throughout 
acquisition  results  in  a  more  spatially  invariant  system  re¬ 
sponse  with  range.  Some  systems  have  Finite  Impulse  Re¬ 
sponse  (FIR)  filters  on  each  receive  channel  to  apply  focal 
time  delay  increments  that  are  smaller  than  the  system 
sampling  interval.  These  filters  usually  have  few  taps  and 
fixed  coefficients.  Among  other  topics,  this  paper  consid¬ 
ers  what  might  be  possible  if  flexible  dynamic  FIR  filters 
were  placed  on  each  channel.  As  we  will  show,  such  an  al¬ 
teration  would  enable  tremendous  control  over  the  system 
response. 

Ultrasound  beam  characteristics  fundamentally  affect 
image  quality  and  the  quality  of  data  acquired  for  signal 
processing.  Because  of  this  significance,  much  of  system 
design  is  dedicated  to  optimizing  beam  parameters.  The 
two  most  important  beam  parameters  are  mainlobe  width 
and  sidelobe  levels.  Mainlobe  width  determines  the  system 
point  resolution,  and  sidelobe  levels  determine  the  system 
contrast.  As  mentioned  previously,  these  and  other  beam 
parameters  can  be  adjusted  by  changing  the  magnitude 
and  phase  (or  time  delay)  of  the  weightings  applied  to 
the  active  elements.  These  parameters  are  also  influenced 
by  the  size  of  the  active  aperture  (the  number  of  active 
elements)  and  the  frequency  of  the  ultrasound  pulse. 
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Magnitude  weightings,  phase  weightings,  aperture  size, 
and  operating  frequency  can  each  be  adjusted  to  manipu¬ 
late  beam  parameters.  Unfortunately,  these  bearaforming 
parameters  do  not  act  independently;  altering  one  changes 
the  impact  of  each  of  the  others.  Thus,  beamforraer  de¬ 
sign  is  a  complicated  multiparameter  optimization  prob¬ 
lem.  Because  of  this  complexity,  beamformer  parameters 
are  typically  determined  using  a  combination  of  ad  hoc 
methods,  simplified  theory,  and  iterative  simulation  and 
experimentation.  Although  these  methods  are  effective, 
they  are  time  consuming  and  provide  no  guarantee  that 
the  optimal  solution  has  been  found.  Therefore,  there  is 
a  fundamental  need  to  develop  a  design  method  that  si¬ 
multaneously  considers  the  impact  of  all  the  controllable 
beamformer  parameters  in  a  straightforward  and  rigorous 
way. 

We  propose  a  general  aperture  design  method,  sup¬ 
ported  by  rigorous  theory,  that  can  be  applied  in  arbitrary 
3ystem  geometries  to  design  apertures  that  optimize  beam 
parameters.  Our  technique  utilizes  a  linear  algebra  formu¬ 
lation  of  the  sum  squared  error  (SSE)  between  the  system 
point  spread  function  (psf)  and  the  desired  or  goal  psf. 
Minimization  of  this  error  yields  unique  aperture  weights 
ings  that  maximize  the  system  psf’s  resemblance  to  the 
desired  psf.  A  strength  of  our  approach  is  that  it  utilizes 
system  characteristics  that  may  be  obtained  through  the¬ 
ory,  simulation,  or  experiments.  Simulation  results  are  pre¬ 
sented  in  another  paper  [2]. 

Our  method  is  similar  to  the  technique  used  by  Ebbini 
and  Cain  [3]  to  generate  specialized  beam  patterns  for  hy¬ 
perthermia,  and  by  Li  et  al  [4]  for  the  compensation  of 
blocked  elements.  However,  there  are  significant  differences 
between  these  methods  and  the  technique  we  describe. 
Chief  among  these  is  the  fact  that  our  technique  is  more 
general.  It  describes  a  method  for  the  design  of  optimized 
apertures  for  any  application,  and  enables  the  design  of 
arbitrary  beam  patterns.  Such  patterns  are,  of  course,  fun¬ 
damentally  constrained  by  wavelength  and  aperture  size. 
The  designed  beam  simply  will  be  that  which  most  closely 
approximates  the  desired  beam  given  the  limitations  of 
the  physics.  The  differences  between  the  MSSE  technique 
and  the  techniques  described  in  [3]  and  [4]  are  discussed 
in  more  detail  later  in  this  manuscript.  This  paper  out¬ 
lines  the  theoretical  description  of  the  MSSE  technique  for 
narrowband  and  broadband  systems,  describes  a  modified 
technique  for  reduced  computational  cost,  and  discusses  a 
few  examples  of  applications.  Simulation  results  for  these 
examples  arc  described  in  [2]. 


II.  Theory 

We  present  one-way  and  two-way  continuous  wave 
(CW)  and  broadband  formulations  of  the  MSSE  technique. 
Please  note  that  all  operators  and  variables  in  the  CW  for¬ 
mulations  are  complex  valued,  and  operators  and  variables 
in  the  broadband  formulations  are  real  valued. 


A.  Continuous  Wave  Formulation 

1.  One-Way  Analysis:  The  phase  and  magnitude  of 
the  ultrasonic  field  at  a  point  in  space  generated  by  an 
ultrasound  transducer  element  depend  upon  several  fac¬ 
tors.  These  include  the  Euclidean  distance  between  the 
point  and  the  element,  the  orientation  of  the  element  rel¬ 
ative  to  the  point,  the  frequency  of  the  emitted  wave,  and 
frequency-dependent  attenuation  of  the  medium.  Assum¬ 
ing  linear  propagation,  the  complex  ultrasonic  field  p(x*) 
at  a  point  in  space  can  be  expressed  as: 

oo 

p(xi)=  J  S(x,Xi)W(x)dX,  (1) 

—  OO 

where  x  represents  position  in  the  aperture  plane,  W(x) 
is  the  complex  aperture  weighting  function,  and  S(x,£,) 
is  a  propagation  function  that  incorporates  any  or  all  of 
the  factors  mentioned  above.  S(x»Zt)  determines  the  com¬ 
plex  field  at  Xi  due  to  the  aperture  weighting  at  x*  This 
propagation  function  may  be  determined  through  theory, 
simulations,  or  experiments.  A  simple  propagation  func¬ 
tion  can  be  formulated  based  on  the  Rayleigh-Sommerfeld 
diffraction  equation,  which  is  derived  in  ([5,  pp.  46-50]). 
The  propagation  function  may  also  include  limited  element 
angular  response  using  the  formulation  derived  in  [6],  as 
well  as  other  such  complicating  factors. 

The  formulation  for  the  ultrasonic  field  can  readily  be 
discretized  because  the  transducer  aperture  comprises  a 
finite  number  of  elements.  Therefore,  the  field  pj  at  a  point 
j  due  to  an  aperture  of  N  elements  can  be  expressed  as: 

N 

Pj  =  J2Si’iWi’  (2) 

1-1 

where  sitj  is  the  propagation  function  that  determines  the 
field  at  the  point  j  due  to  the  ith  element,  and  is  the 
weighting  applied  to  the  ith  element.  Therefore,  pj  is  the 
sum  of  the  contribution  of  each  element  to  the  field  at 
j.  Note  that  this  formalism  makes  no  assumptions  about 
the  geometry  or  other  characteristics  of  the  transducer  el¬ 
ements.  Eq.  (2)  can  be  expressed  in  a  matrix  formulation 
as  follows: 


Pj  —  [sij  82, j  Sz,j  *  ' 


W  i 
W2 
Wz 


U)  AT 


(3) 


Therefore,  the  one-way  Af-point  lateral  psf  at  the  range 
z  can  be  represented  as: 


P*  = 


Sl,l  $2,1  •••  S/v,  1 

$1, 2  S2>2  **' 

w2 

*  **• 

WN 

SZW,  (4) 
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where  Sz  is  an  Af  x  N  matrix  of  propagation  functions  with 
Sij  denoting  the  propagation  function  that  determines  the 
field  at  the  point  j  due  to  the  ith  element,  W  is  an  N  x  1 
vector  of  aperture  weightings,  and  the  resulting  psf  Pz  is 
an  M  x  1  vector.  Note  that  our  field  point/element  no¬ 
tation  differs  from  the  conventional  row/column  ordering. 
We  have  chosen  this  ordering  for  purposes  of  clarity  and 
consistency  within  the  manuscript,  as  will  become  more 
clear  in  the  broadband  formulation  below.  As  mentioned 
previously,  this  formulation  permits  analysis  with  compli¬ 
cated  propagation  functions  that  may  include  limited  ele¬ 
ment  angular  response,  frequency  dependent  attenuation, 
and  other  factors  that  are  difficult  to  model.  Using  the  for¬ 
mulation  in  (4),  the  transmit  and  receive  psfs  at  the  range 
z  can  be  expressed  as: 


important  to  enforce  low  sidelobe  levels  than  to  precisely 
control  the  mainlobe.  In  such  cases,  we  can  incorporate  a 
weighting  function,  F ,  that  emphasizes  or  de-emphasizes 
selected  features  in  the  psf  during  the  MSSE  design  pro¬ 
cess.  The  SSE  (7)  can  be  rewritten  with  the  weighting 
function  as: 

SSE  =  (Fd SZT  —  Fd  Prx) H  (Fd SZT-Fd  PTz)  , 

where  Fd  is  a  diagonalized  MxM  matrix  with  the  elements 
of  F  along  its  0th  diagonal.  These  elements  should  have 
a  large  value  where  a  close  match  between  the  goal  and 
designed  psf  is  required.  Smaller  values  can  be  .used  in 
regions  where  a  close  match  is  less  critical.  The  solution 
for  the  receive  weightings,  as  drawn  from  [7]  is: 


Pt*  =  SzT  (5) 

and 

Frz  ==  R,  (6) 

respectively,  where  T  and  R  are  the  transmit  and  receive 
aperture  weightings,  respectively. 

Let  Ftz  represent  the  desired  one-way  psf  for  the  ap¬ 
plication  of  interest.  We  can  then  characterize  the  degree 
of  similarity  of  the  desired  psf,  Ptz ,  and  the  actual  sys¬ 
tem  psf,  Ptz,  by  the  sum  squared  error  between  them. 
Minimizing  this  SSE  would  yield  a  system  psf  optimally 
similar  to  the  goal  psf.  Therefore,  beamformer  design  is 
simply  the  selection  of  transmit  aperture  weightings  such 
that  the  SSE  between  the  desired  and  actual  system  psfs  is 
minimized.  Using  (5),  the  SSE  can  be  expressed  as  follows: 


T  =  ((FdS,f  FdSz)  1  (FdSzf  Fd  PTz 

=  (Sf  Fd  FdSz) -1  S»  Fd  Fd  Prz  {1Q} 

=  (FdSzf  FdPrz, 

where  (Fd  Sx)#  is  the  pseudoinverse  of  FdSz. 

$.  Two-  Way  Analysis:  In  most  ultrasonic  imaging  ap¬ 
plications,  the  two-way  impulse  response  is  of  greater  inter¬ 
est  than  the  one-way  response.  The  two-way  response  can 
be  readily  determined  by  applying  the  radar  equation  [8]. 
It  states  that,  for  continuous  wave  applications,  the  two- 
way  response  is  simply  the  product  of  the  transmit  and 
receive  responses.  Applying  this  knowledge  to  our  linear 
algebra  formulation  yields  the  following  two-way  response: 


SSE  =  (Ptz  -  Ptz')  {Ftz  -  Ptz) 

—  (sz  T  -  Prz)H  (sz  T  -  PTz )  , 


Ptrz  =  Ptz  ■  Prz  =  (S,  T)  •  ( Sz  R) ,  (11) 

where  indicates  point-by-point  multiplication.  Eq.  (11) 
can  be  rewritten  as: 


where  the  superscript  H  denotes  a  conjugate  transpose 
operation. 

The  formulation  in  (7)  is  common  in  signal  process¬ 
ing,  and  substantial  literature  is  devoted  to  the  solution  to 
(7)  with  the  minimum  SSE  (least  squares  solution).  Draw¬ 
ing  upon  [7],  the  optimal  transmit  aperture  weightings  are 
given  by: 

T  =  {Si1  Szyl  Sf  Ptz  =  sf  Ptz,  (8) 

where  the  superscripts  ‘—l*  and  denote  a  matrix  in¬ 
verse  and  a  pseudoinverse  operation,  respectively.  There¬ 
fore,  Sf  is  the  pseudoinverse  of  Sz. 

Eq.  (8)  provides  a  simple  method  for  the  calculation  of 
the  transmit  weightings  that  yield  the  system  psf  at  the 
range  z  that  is  optimally  similar  to  the  goal  or  desired  psf. 

2.  One-  Way  Analysis  with  Weighting  Function:  In  cer¬ 
tain  applications,  the  psf  characteristics  at  some  lateral 
positions  may  be  more  critical  than  at  others  because  the 
width  of  the  mainlobe  of  the  psf  determines  system  point 
resolution  and  the  magnitudes  of  the  side-lobes  determine 
the  system  contrast.  In  a  given  application,  it  may  be  more 


Ptrz  =  Przd  Prz  =  Ptza  SzR  —  Ptzas  R> 

(12) 


where  PrZd  is  a  diagonal  MxM  matrix  with  the  elements 
of  Ptz  along  its  0th  diagonal,  and  Przds  =  PrzdSz-  This 
changes  the  point  multiplication  operation  to  a  regular 
matrix  multiplication  operation. 

Similar  to  the  one-way  analysis,  if  the  goal  two-way  psf 
for  the  application  of  interest  is  Ptrz,  the  SSE  can  be 
expressed  as: 

(13) 


SSE  =  ^ Ptrz  -  Pr Hz')  (Ptrz  —  Ptrz ) 


Applying  (12),  (13)  can  be  rewritten  and  solved  using  [7] 
as  shown  below,  yielding  the  optimum  receive  weightings: 

SSE  =  (Przds  R  -  Ptrz^  (Przds  R  -  Ptrz )  , 

R  ~  ( pTzdS  PTzdS )  PrzdS  Ptrz  =  P *TzdS  PtRz > 

(15) 


where  Pfzds  is  the  pseudoinverse  of  Przds • 
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Eq.  (15)  specifies  the  complex  weightings  to  be  applied 
to  the  transducer  elements  constituting  the  receive  aper¬ 
ture  to  obtain  a  system  psf,  PTrz ,  at  the  range  2,  that 
optimally  resembles  the  desired  or  goal  psf,  Ptrz ■ 

4 •  Two-Way  Analysis  with  Weighting  Function:  As 
with  the  one-way  analysis,  a  weighting  function  F  can  be 
incorporated  in  the  SSE  formulation  to  selectively  empha¬ 
size  or  de-emphasize  features  in  the  psf  during  the  mini¬ 
mization  operation.  Rewriting  the  SSE  (14)  after  including 
the  weightings  matrix,  the  resulting  receive  weightings  can 
be  solved  for  in  a  manner  similar  to  that  used  in  (9)  and 
(10). 

R  =  ((F<f  Przd$)H  Fd  Przdsj  ( Fd  Przds )H  Pd  Ptrz 

=  {P&is  Fd  Fd  PTzds)~l  p?zds  F%  Fd  PTRt 
=  ( Fd  Przds )*  FdPrRz, 

where  Fd  is  again  a  diagonalized  M  x  M  matrix  with  the 
elements  of  F  along  its  0th  diagonal,  and  {Fd  Przds)#  is 
the  pseudoinverse  of  Fd  Przds  • 

B.  Broadband  Formulation 


1.  One-Way  Analysis:  The  formulations  described  in 
the  previous  section  are  CW  formulations.  Medical  ultra¬ 
sound  systems  use  monochromatic  (CW)  excitation  only 
for  specific  modalities  such  as  CW  Doppler.  For  the  ma¬ 
jority  of  applications,  an  ultrasound  pulse  with  some  finite 
bandwidth  is  used.  The  previously  described  CW  formula¬ 
tion  will  have  limited  accuracy  in  this  broadband  scenario. 

The  one-way  psf  at  a  specific  range  z  can  be  represented 
as  a  function  of  time  and  lateral  position  as  follows: 


Pl.l  P2,l 
Pi, 2  p2,2 


|Pl,ntp  * 


Pnp 


(17) 


where  Pz  is  an  ntp  x  np  matrix  that  is  a  two-dimensional 
function  of  position  and  time.  It  consists  of  the  field  at 
each  of  Up  lateral  points  in  space,  at  each  of  ntp  points  in 
time.  It  comprises  elements  of  the  form  pitjt  which  is  the 
field  at  lateral  point  i  at  time  sample  j.  Eq.  (17)  can  be 
•  rewritten  by  reshaping  the  matrix  as  follows: 


Pi,i 
Pi, 2 


Pl,ntp 

P2,l 

P2,2 


[Pnp,ntpJ 


(18) 


where  the  first  ntp  elements  represent  the  field  at  lateral 
point  1  at  each'  of  nip  time  samples,  the  next  7itp  elements 
represent  the  field  at  lateral  point  2  for  the  same  ntp  time 
samples,  and  so  on  until  the  last  ntp  time  samples  for  the 
riph  lateral  point. 

The  field  at  a  point  in  space  over  ntp  time  samples  can 
be  expressed  as  a  function  of  a  propagation  matrix,  A,  and 
a  set  of  aperture  weightings,  T.  The  propagation  matrix 
depends  upon  the  excitation  pulse  and  the  element  impulse 
responses  of  the  transmit  aperture.  It  is  a  function  of  time 
and  the  spatial  positions  of  the  element  and  field  point 
under  consideration.  It  describes  the  contribution  of  each 
element  at  each  field  point  as  a  function  of  time.  The  aper¬ 
ture  weightings  are  also  two-dimensional,  being  a  function 
of  the  element  number  and  time,  and  can  be  expressed  for 
each  of  na  elements  over  each  of  nta  time  samples  as: 


T  = 


*1,1  $2,1 
$1,2  *2,2 


1*1*"*® 


t 


no  ,nJa 


(19) 


where  Uj  is  the  aperture  weighting  for  element  i  at  time 
j •  These  weightings  essentially  form  the  coefficients  of  a 
FIR  filter.  Eq.  (19)  also  can  be  reshaped  as  follows: 


*1,1 

*1,2 


*1  ,n«a 
*2,1 
*2,2 


(20) 


t 


where  the  first  Uta  elements  axe  the  aperture  weightings 
for  element  1  at  each  of  nta  time  samples,  the  next  nta 
elements  are  the  weightings  for  element  2  for  the  same  nta 
time  samples,  and  so  on  untfi  the  last  nta  elements  for  the 
n‘h  element.  Using  (18)  and  (20),  we  can  now  write  the 
complete  one-way  system  psf  Pz  as  (21)  (see  next  page)  or 


Pz  =  AzT ,  (22) 

where  Az  is  an  (np  •  ntp)  x  (na  •  nta)  propagation  matrix. 
Each  element  ai^i  determines  the  field  due  to  the  weight¬ 
ing  at  time  sample  j  applied  to  element  z,  at  time  sample 
l  at  lateral  field  point  k .  The  transmit  and  receive  psfs  at 
the  range  z,  therefore,  can  be  expressed  as  follows: 


Ptz  =  AzT1  (23) 

and 


Prz  =  Az  R,  (24) 

respectively,  where  T  and  R  are  the  (na  -nta)  x  1  transmit 
and  receive  aperture  weightings,  respectively. 
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Fig.  1.  Broadband  formulation.  The  pulse  applied  to  each  element 
is  convolved  with  a  set  of  weights  that  is  distinct  for  each  element. 
This  operation  is  analogous  to  applying  an  FIR  filter  to  each  ele* 
ment/channel. 


We  can  now  derive  the  transmit  weightings  that  force 
the  one-way  system  psf  to  optimally  resemble  a  specified 
goal  psf,  Ptz-  The  SSE  between  the  goal  and  system  one¬ 
way  psfs  at  range  z  can  be  expressed  as: 

SSE  =  (pTz  -  Ptz)T  (Ptz  -Ptz),  (25) 

where  the  superscript  T  indicates  a  transpose  operation. 
A  conjugate  operation  is  unnecessary  because  the  weights 
are  real  valued.  Substituting  (23)  in  (25)  yields: 

SSE  =  (azT- PTiy  (az T -  PTz)  ■  (26) 

As  in  the  CW  formulations,  the  transmit  weightings 
that  minimize  the  SSE  can  be  determined  using  [7]  and 
are  given  by: 

T  =  (Af  A.)'1  A f  PTz  =  A#  Ptz,  (27) 
where  Af  is  the  pseudoinverse  of  Az. 

2.  One-  Way  Analysis  with  Weighting  Function :  As  in 
the  CW  analysis,  we  can  incorporate  a  weighting  function 
F  in  the  analysis.  In  the  broadband  case,  however,  F  is  an 
(ftp  ntp)  xl  element  weighting  vector  as  shown  below. 

f  1 1 


The  generation  of  a  one-way  psf  is  shown  in  Fig.  1.  The 
transmit  pulse  applied  to  each  element  is  first  convolved 
with  a  set  of  weights  that  is  unique  for  each  element.  The 
results  of  these  convolutions  then  are  convolved  with  the 
respective  element  spatial  impulse  responses  and  finally 
summed  in  space  to  generate  the  beam  pattern.  Note  that 
the  broadband  formulation  in  (22)  describes  the  psf  at  a 
single  range.  Therefore,  the  formulation  does  not  perform 
a  convolution  to  determine  the  ultrasonic  field  at  multi¬ 
ple  range  points.  However,  there  is  a  convolution  that  de¬ 
scribes  the  response  at  each  lateral  point  of  interest  as  a 
function  of  time  in  the  formulation.  The  transmit  weights 
that  are  calculated  using  the  formulation  are  eventually 
convolved  with  the  pulse  (and  the  element  response)  to 
form  the  resultant  beam  pattern.  This  operation  would 
impact  the  field  at  multiple  ranges. 


hi  ’ 

/2,2 


(28) 


\Jnv,ntp 


where  F  consists  of  the  weighting  to  be  applied  to  the  field 
at  each  of  np  points  in  space,  at  each  of  ntp  points  in  time. 
It  comprises  elements  of  the  form  faj,  which  is  the  weight¬ 
ing  applied  to  the  field  at  point  t  at  time  sample  j.  The 
SSE  (26)  can  be  rewritten  with  the  weighting  function  as: 


SSE  =  (Fd  AZT  -  Fd Pt*)T  (E,  AZT  —Fd PTz)  , 


(29) 
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where  Fd  is  a  diagonalized  (np  •  ritv)  x  (np  •  rifP)  matrix 
with  the  elements  of  F  along  its  0th  diagonal.  The  solution 
for  the  transmit  weightings,  therefore,  is: 

T=({FdAz)T  FdAz)~l  (FdAzf  FdPTz 

~  Fd  At)  1  A?  Fd  Fd  PTz 
=  {FdAzfFdPTz ,  (30) 

where  (FdAz)#  is  the  pseudoinverse  of  FdAz. 

3.  Two-Way  Analysis :  The  two-way  pulse-echo  psf  can 
also  be  expressed  as  a  linear  algebra  formulation  in  a  sim¬ 
ilar  fashion  to  the  one-way  psf: 


Ptrz  —  Azz  Ry  (31) 

where  Azz  is  the  propagation  function,  and  R  is  the 
(fta  *  n*a)  *  1  weighting  vector  with  weights  for  each  of 
na  receive  elements  at  each  of  nta  time  samples.  Azz  is  a 
function  of  the  transmit  aperture  weights,  the  excitation 
pulse,  and  the  element  impulse  responses  of  the  transmit 
and  receive  apertures.  The  SSE  between  the  goal  and  sys¬ 
tem  pulse-echo  psfs  at  range  z  is: 

SSE  =  (Ptrz  -  Pti tr)  (Pthz  -  PtrJ)  , 

'  (32) 

where  Ptrz  is  the  goal  or  desired  pulse-echo  psf  at  range 

z.  We  can  substitute  (31)  in  (32)  and  solve  for  the  receive 
weightings  to  be  applied  using  (7]: 

SSE  =  (azz  R  -  Ptkz)  (azz  R  -  Ptrz)  , 

'  (33) 

P  ~  {A^z  Azz )  A^z  Ptrz  =  Afz  PTrz , 

(34) 

where  Afz  is  the  pseudoinverse  of  Azz. 

4-  Two-Way  Analysis  with  Weighting  Function :  We 
,  again  can  include  a  weighting  function  F,  which  is  an 
{rip  •  ntp)  x  1  element  weighting  vector.  Solving  after  appli- 
,  cation  of  the  weighting  function  yields  the  following  receive 
weightings: 

R—  FdAzz )  FdAzz'j  (Frf  AZZ)T  PdPTRz 

=  {ATZZ  Fj Fd  Azz)'1  AJZ  Fj  Fd  PTRz 
=  (FdAzzf  FdPTRzf  (  j 


where  Fd  is  again  a  diagonalized  (np  •  ntp)  x  (np  *  ntp)  ma¬ 
trix  with  the  elements  of  F  along  its  0th  diagonal,  and 
{Pd  Azz )#  is  the  pseudoinverse  of  Fd  Azz. 


technique  requires  significant  resources,  due  to  the  pseu¬ 
doinverse  operation  and  the  large  propagation  matrices  re¬ 
quired.  Note  that  the  application  of  the  weights  has  a  much 
lower  computational  cost  than  the  design  of  the  weights. 
In  order  to  reduce  the  computational  complexity  of  our 
broadband  formulation,  we  take  advantage  of  the  symme¬ 
try  present  in  the  system.  We  first  use  the  lateral  symme¬ 
try  of  the  psf.  This  symmetry  means  that  we  can,  if  we 
choose,  use  just  half  of  both  the  goal  and  system  psfs  for 
the  calculation  of  the  optimal  weightings.  The  one-way  psf 
is  then: 


P\,l  P2,l 

'  Pn p/2,1 

Pi, 2  P 2,2 

:  : 

Pi, ntp  • 

Pnp/2,fitp_ 

where  Pz  is  an  (np/2)  x  (ntp)  matrix,  consisting  of  the  field 
at  each  of  only  rip/2  points  in  space  at  each  of  ntp  points 
in  time.  Eq.  (36)  can  be  rewritten  as: 


Pi,i 
Pi  ,2 


ft  = 


Pl,ntp 

P2,l 

P2,2 


LPnP/2,ntp  J 


(37) 


where  Pz  is  an  (np/2  •  ntp)  x  1  vector.  The  goal  psf  also 
must  be  rewritten  in  a  similar  fashion. 

The  symmetry  of  the  transmit  and  receive  apertures  is 
another  property  that  can  be  exploited  to  reduce  compu¬ 
tational  cost.  As  shown  in  Fig.  2,  pairs  of  elements  can  be 
generated  by  grouping  elements  on  either  side  that  are  at 
the  same  distance  from  the  center  axis,  because  these  ele¬ 
ments  will  have  the  same  weightings,  assuming  no  beam¬ 
steering.  Assuming  that  the  aperture  comprises  an  even 
number  of  elements,  the  transmit  weightings  then  can  be 
expressed  as: 


*1,^,1 

^l,n0>2 

<2,na-l,l  * 
*2,na~l,2  * 

*  *n0/2,n0/2+l,l 

t\ tn«,nla 

* 

‘  ^nQ/2,nQ/2-fltnta. 

(38) 


5.  Reduced  Computational  Cost  Through  Symmetry  Re¬ 
lations:  The  computation  of  aperture  weights  in  the  MSSE 


where  is  the  aperture  weighting  for  elements  i  and  j 
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Fig.  2.  Exploitation  of  symmetry  for  reduced  computational  cost. 
Because  the  aperture  is  symmetric  about  the  center  elements,  pairs 
of  elements  that  would  have  the  same  weights  are  grouped  together. 
AlsOj  the  lateral  symmetry  of  the  psf  about  the  center  axis  permits 
analysis  with  just  one  half  of  the  psf. 


at  time  fc.  Eq.  (38)  also  can  be  reshaped  as  follows: 


£l,na,l 

^l,na»2 


T  = 


tl,n0,nt0 

$2,na-l,2 


L^na/2,na/2+l,n«aJ 


(39) 


where  the  first  nta  elements  are  the  aperture  weightings 
fox  elements  1  and  na  at  each  of  nta  time  samples,  the 
next  nta  elements  are  the  weightings  for  elements  2  and 
na  —  1  at  the  same  nta  time  samples,  and  so  on  until  the 
last  nta  weightings  for  elements  na/2  and  naj 2  -  1. 

The  derivation  is  now  analogous  to  the  analysis  in  (23) 
and  (25)-(27)  yielding  the  transmit  weightings: 

T  =  (Atz  A*) -1  Atz  Pt*  =  Af  PTz.  (40) 


This  use  of  symmetry  reduces  the  size  of  the  matrix 
subject  to  the  pseudoinverse  operation  by  a  factor  of  2 
in  each  dimension;  therefore,  it  also  reduces  the  required 
memory.  The  same  technique  can  be  applied  in  the  two- 
way  pulse-echo  scenario  to  reduce  computational  cost. 


III.  Applications 

The  MSSE  technique  described  above  is  extremely  gen¬ 
eral  and  can  be  applied  in  wide-ranging  scenarios.  A  few 
possible  applications  are  described. 

A.  Enhanced  Depth  of  Field 

The  DOF  of  an  ultrasound  imaging  system  is  generally 
defined  as  the  axial  region  over  which  the  system  is  in  fo¬ 
cus,  or  more  rigorously,  the  axial  region  over  which  the 
system  response  satisfies  some  chosen  criterion.  It  is  gen¬ 
erally  desired  that  the  system  psf  remains  similar  to  the 
psf  at  the  focus  for  as  large  an  axial  span  as  possible. 

Current  techniques  to  improve  DOF  include  transmit 
apodization,  dynamic  receive  apodization,  and  dynamic 
receive  focusing  (1].  However  effective  the  above  techniques 
are,  their  implementation  is  typically  ad  hoc  and  lacks 
formal  theory  describing  their  effectiveness  in  improving 
DOF.  If  the  MSSE  technique  is  implemented  for  every 
range  under  consideration,  with  the  goal  psf  being  the  psf 
obtained  at  the  focus,  we  can  derive  formally  the  receive 
apodization  weightings  that  force  the  psf  at  each  specific 
range  of  interrogation  to  be  maximally  similar  by  min¬ 
imizing  the  SSE.  These  weightings  can  then  be  used  to 
implement  dynamic  apodization  and  maximize  the  DOF. 

As  demonstrated  in  (11)  and  (12),  we  can  express  the 
pulse-echo  two-way  psf  at  a  range  z  as  the  point- by-point 
multiplication  of  the  one-way  transmit  and  receive  psfs  at 
range  z\ 

Ptrs  =  Pt z  •  Prz  =  (S.  T)  -  (SM  R)  '  (41) 

PtRz  =  PTzd  PRz  “  Przd  Sz  R  =  PrzdS  R • 

(42) 

The  transmit  psf,  Ptz,  is  fixed  at  each  range  of  interro¬ 
gation  because  we  consider  each  transmit  focus  separately. 
From  (15),  the  receive  weightings  at  a  range  z  that  mini¬ 
mize  the  SSE  between  the  psf  at  the  focus  and  the  psf  at 
z,  and  therefore  maximize  the  DOF  are: 

R  —  (P-TidS  PTzds)  Ptz&S  PtRz  =  P^zdS  PtRzi 

(43) 


where  Ptrz  is  the  psf  at  the  focus. 

B.  Limited  Diffraction  Transmit  Beams 

Modem  ultrasound  systems  use  dynamic  receive  focus¬ 
ing  and  dynamic  apodization  to  expand  their  useful  DOF. 
Ultimately,  however,  the  DOF  is  limited  by  the  use  of  a 
fixed  transmit  focus.  Thus,  multiple  transmissions  with 
different  focal  ranges  must  be  performed  along  each  image 
line  in  order  to  obtain  a  high  quality  image.  This  slows 
image  acquisition  significantly. 

Limited  diffraction  beams  [9]  have  been  suggested  as  a 
way  to  enhance  DOF  without  requiring  multiple  transmis¬ 
sions.  The  MSSE  technique  described  here  can  be  readily 
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applied  to  design  such  beams.  We  present  CW  analysis  to 
generate  limited  diffraction  beams  through  the  application 
of  appropriate  transmit  weights. 

If  Pi,p2,P3i  •••  j Pq  are  CW  transmit  (one-way)  psfe 
at  Q  different  ranges  of  interest,  and  we  axe  interested 
in  maintaining  the  transmit  beam  profile  through  these 
ranges,  the  objective  is  to  derive  a  set  of  transmit  aper¬ 
ture  weightings  that  would  minimize  the  SSE  between  the 
system  and  goal  psfe  at  each  of  these  ranges.  We  can  ex¬ 
press  the  one-way  psfe  at  all  ranges  of  interest  as  follows: 


'/V 

fSil 

p2 

s2 

Pz 

—  Ptq  ~ 

S3 

Pq. 

Sq. 

m =sqt, 


(44) 


where  Ptq  is  a  (M  *  Q)  x  1  vector  made  up  of  Q  vertically 
tiled  M-point  psfe  at  the  Q  ranges  of  interest,  Sq  is  an 
(Af  •  Q)  x  N  matrix  made  up  of  Q  vertically  tiled  M  x 
N  propagation  matrices,  one  for  each  of  the  Q  ranges  of 
interest,  and  T  is  the  N  x  1  vector  of  transmit  aperture 
weightings.  The  SSE  then  can  be  expressed  as  follows: 


SSE  —  -  Ptq}  (Ptq  -  Ptq )  , 


where  Ptq  — 


Pim 

P2 

h 


Uw 


(45) 


is  an  (Af  Q)  x  1  vector  consisting  of  Q  vertically  tiled  Af  x  1 
goal  or  desired  psfe,  one  for  each  range.  Note  that  Ptq 
cannot  be  constructed  by  simply  replicating  the  same  psf 
Q  times,  because  constructing  one-way  psfe  with  the  same 
phase  at  each  range  is  impossible.  In  CW  implementation, 
an  appropriate  phase  term  that  accounts  for  propagation 
must  be  applied  to  the  goal  psf  at  each  range.  This  phase 
term  is  a  function  of  the  range  under  consideration  and  the 
wave  propagation  speed.  In  broadband  implementation, 
appropriate  time  delays  need  to  be  applied.  Using  (44), 

;  (45)  can  be  rewritten  and  solved  to  obtain  the  transmit 
weightings  that  minimize  the  SSE: 

SSE  =  (SqT  —  Ptq ^  (^SqT  —  Ptq ^  ,  (46) 

P  *=  (Sq  Sq)  Sq  Ptq  =  Sq  Ptq ,  (47) 

*  where  Sq  is  the  pseudoinverse  of  Sq . 

Eq.  (47)  specifies  the  complex  weightings  to  be  applied 
to  the  transmit  aperture  to  obtain  one-way  transmit  psfe 
at  the  specified  ranges  that  optimally  resemble  the  goal 
.  psfe.  Therefore,  these  transmit  weightings  will  result  in  a 
limited  diffraction  transmit  beam,  which  will  have  a  ben¬ 
eficial  impact  on  DOF. 


Fig.  3.  Reduced  DOF  in  translated  aperture  geometries,  (a)  depicts 
the  angular  scatter  geometry,  and  (b)  shows  the  backscatter  geome¬ 
try.  The  reduced  depth  of  field  is  due  to  the  limited  region  of  overlap 
of  the  transmit  and  receive  beams  in  (a),  as  compared  to  the  com¬ 
pletely  coincident  transmit  and  receive  beams  in  (b). 


C .  Increased  Correlation  Depth  of  Field 
in  Translated  Aperture  Geometries 

Biological  tissues  are  known  to  exhibit  variations  in  an¬ 
gular  scattering  [10]— [12].  That  is,  the  scattered  echo  mag¬ 
nitude  and  phase  depend  upon  the  angle  between  the  prop¬ 
agation  vectors  of  the  incident  and  observing  waves.  It  has 
long  been  hypothesized  that  this  parameter  might  offer 
valuable  diagnostic  information.  We  have  proposed  using 
the  translating  apertures  algorithm  (TAA),  as  the  foun¬ 
dation  of  angular  scatter  imaging  methods  [13].  Previous 
methods  that  were  used  to  make  angular  scatter  measure¬ 
ments  [14],  [15]  have  entailed  the  use  of  pistons  that  were 
mechanically  rotated  to  measure  the  average  angular  scat¬ 
ter  over  some  area  at  a  single  frequency.  Imaging  systems 
also  have  been  developed  to  make  images  at  a  single  scat¬ 
tering  angle  other  than  180°  [16],  [17].  However,  in  all  these 
methods  the  speckle  pattern  that  was  obtained  varied  with 
angle,  and  images  obtained  at  different  angles  could  not 
be  processed  coherently  to  obtain  accurate  complex  angu¬ 
lar  scatter  information.  The  reason  for  the  change  in  the 
speckle  pattern  was  a  rotation  of  the  system  psf  with  a 
change  in  the  relative  position  of  the  transmit  and  receive 
apertures.  In  the  TAA,  the  transmit  and  receive  apertures 
are  translated  by  an  equal  distance  in  opposite  directions. 
This  enables  the  acquisition  of  accurate  angular  scatter 
data  without  the  confounding  influence  of  system  psf  ro¬ 
tation. 

Although  offering  several  advantages,  the  TAA  results 
in  a  significantly  reduced  DOF  as  the  transmit  and  receive 
apertures  are  translated.  This  is  due  to  the  crossing  of  the 
transmit  and  receive  beams,  which  results  in  interference 
over  a  reduced  area  as  the  apertures  are  translated.  This 
is  shown  in  Fig.  3.  Our  technique  of  optimizing  dynamic 
receive  aperture  weightings  can  be  applied  to  improve  the 
correlation  DOF  between  the  backscatter  (nontranslated) 
and  angular  scatter  (translated)  geometry  psfe  at  a  specific 
range.  Note  that  the  previously  described  symmetry  tech¬ 
nique  for  the  broadband  formulation  cannot  be  completely 
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utilized  in  translated  aperture  geometries.  This  is  due  to 
a  loss  of  symmetry  that  is  caused  by  the  translation  of  the 
apertures,  which  results  in  unique  element  weightings  for 
each  element.  We  can  write  the  two-way  CW  psf  for  the 
backseat  ter  geometry  at  range  z  as  follows: 

Ptrzo  —  Ptzo  •  PrzO  =  (Stzo  To)  •  {SrzO  Rq)  } 

(48) 

where  Ptzo  and  PrzO  are  the  M  x  1  transmit  and  receive 
psfs,  respectively;  SrzO  and  SrzO  are  the  M  x  N  transmit 
and  receive  propagation  functions,  respectively,  at  range 
z;  and  To  and  Rq  are  the  N  x  1  transmit  and  receive  aper¬ 
ture  weightings,  respectively.  The  subscript  ‘O’  denotes  no 
translation  (zero  shift),  or  the  backscatter  geometry.  Simi¬ 
larly,  the  two-way  CW  psf  for  the  angular  scatter  geometry 
at  the  same  range  z  is: 

Ptrzi  =  Ptz\  *  Przi  —  ( Stz\  Ti)  *  (SW  Ri) , 

(49) 

where  Ptz i  and  Prz\  are  the  M  x  1  transmit  and  receive 
psfs,  respectively;  Stzi  and  Srz i  are  the  M  x  N  transmit 
and  receive  propagation  functions,  respectively,  at  range 
z;  Ti  and  R\  are  the  N  x  1  transmit  and  receive  aperture 
weightings,  respectively;  and  the  subscript  ‘1’  denotes  the 
translated  angular  scatter  geometry.  Because  the  apertures 
are  translated,  the  propagation  functions  are  no  longer  the 
same  for  the  transmit  and  receive  apertures.  Eq.  (49)  can 
be  rewritten  as: 

PtRzI  =  P'fzld  PrzI  =  Przld  SrzI  R\  =  PrzldS  Rli 

(50) 

where  PrzU  is  a  diagonalized  M  x  M  matrix  with  the 
elements  of  Przi  along  its  0th  diagonal,  and 

PtzUS  =  Przid  Srzi.  (51) 

Our  objective  is  to  maintain  a  constant  system  response 
as  the  apertures  are  translated.  Again,  this  means  maxi¬ 
mizing  the  correlation  between  the  system  responses  of  the 
backscatter  and  angular  scatter  geometries,  or  minimizing 
the  SSE  between  the  two  psfs.  The  SSE  can  be  expressed 
as  follows: 

SSE  =  ( Ptrzi  ~  Ptrzo)H  ( Ptrzi  —  Ptrzo)  • 

(52) 

Using  (50),  (52)  can  be  modified  and  solved  to  obtain 
the  receive  weightings  that  minimize  the  SSE  as: 

SSE  =  ( Przids  Ri  -  Ptrzo)H  ( Przids  Ri  —  Ptrzo)  » 

(53) 

Rl  =  ( PrzldS  Przids)  P?zldS  PrRzO  =  PrzldS  PtRzO* 

(54) 

where  Pfz\ds  *s  pseudoinverse  of  Przids • 

Eq.  (54)  specifies  the  complex  weightings  to  be  applied 
to  the  transducer  elements  that  comprise  the  receive  aper¬ 
ture  after  translating  the  apertures.  Application  of  these 


weightings  would  generate  a  system  psf  Ptrzi  at  the  range 
z  that  optimally  resembles  the  goal  psf  Ptrzo,  or  the  psf 
with  no  translation. 

D.  Optimal  Receive  Weighting  for  Harmonic  Imaging 

The  previously  described  MSSE  technique  is  not  lim¬ 
ited  to  aperture  design  assuming  linear  propagation.  Con¬ 
ventional  ultrasound  imaging  assumes  linear  propagation 
of  the  ultrasound  pulse,  and  the  receive  process  assumes 
that  the  received  echoes  have  the  same  frequency  content 
as  that  of  the  transmitted  pulse.  However,  the  propaga¬ 
tion  process  is  nonlinear  and  shifts  some  of  the  signal  en¬ 
ergy  to  harmonics  of  the  fundamental  frequency.  Current 
state-of-the-art  ultrasound  systems  have  the  capability  of 
imaging  echoes  received  at  these  higher  harmonics  (har¬ 
monic  imaging),  for  improved  image  contrast  and  reso¬ 
lution  [18].  Our  technique  of  dynamic  weighting  can  be 
adapted  to  design  receive  apertures  for  harmonic  imaging. 
The  transmit  beam  profile  resulting  from  nonlinear  pro¬ 
cesses  can  be  determined  analytically,  experimentally,  or 
through  simulations,  and  be  substituted  into  our  formu¬ 
lation.  The  algorithm  assumes  linear  propagation  during 
the  receive  process.  Eq.  (15),  which  is  rewritten  below,  de¬ 
scribes  the  relationship  between  the  optimum  weightings 
and  the  analytically  or  experimentally  determined  goal  psf. 

P  =  ( P?zdS  PrzdS )  P-TzdS  PrRz  =  P^zdS  PtRz  • 

IV.  Discussion 

The  MSSE  technique  is  a  generalized  technique  for  the 
design  of  arbitrary  system  responses,  through  the  use  of 
aperture  weights.  Analysis  that  illustrates  the  theory  un¬ 
derlying  the  technique  has  been  presented  for  both  CW 
and  broadband  systems.  The  propagation  functions  that 
have  been  described  can  be  determined  through  experi¬ 
ments  and  simulations  or  derived  from  theory. 

The  real-time  implementation  of  the  beamformer  is  con¬ 
ceptually  simple.  In  the  CW  case,  complex  apodization 
would  be  used  to  implement  the  technique.  The  weights 
either  can  be  stored  and  retrieved  in  a  look-up  table  or 
calculated  as  required.  For  the  broadband  case,  implemen¬ 
tation  is  analogous  to  using  an  FIR  filter  on  each  channel. 
In  the  one-way  design  implementation,  conventional  FIR 
filters  would  be  used;  in  the  two-way  design  implemen¬ 
tation,  dynamic  shift-variant  FIR  filters  would  be  used. 
There  are  no  restrictions  on  filter  length;  the  algorithm 
uses  the  specified  length  to  optimize  performance.  How¬ 
ever,  as  in  conventional  FIR  filters,  performance  depends 
upon  the  filter  length,  and  it  is  advisable  to  have  the  maxi¬ 
mum  possible  length  that  can  be  practically  implemented. 
Some  current  ultrasound  systems  already  have  a  crude  FIR 
filter  on  each  channel,  and  they  conceivably  could  be  ex¬ 
tended  to  apply  the  MSSE  technique  by  modifying  them 
to  be  shift  variant.  Again,  the  weights  can  either  be  stored 
and  retrieved  in  a  look-up  table  system  or  calculated  as 
required. 
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The  design  method  described  here  is  different  and  su¬ 
perior  to  the  pseudoinverse  based  methods  described  in  [3] 
and  [4]  for  several  reasons.  The  methods  in  [3]  and  [4]  are 
limited  to  CW  systems.  We  present  analysis  for  CW  as  well 
as  broadband  systems.  Our  method  is  also  more  general  as 
we  describe  a  general  method  for  the  design  of  apertures 
for  any  application.  Another  important  distinction  is  that 
the  methods  in  [3]  and  [4]  are  limited  to  one-way  analysis 
(i.e.,  either  transmit  or  receive).  Our  technique  is  adapt¬ 
able  for  one-way  or  two-way  analysis.  Another  distinction 
is  that  in  [3]  and  [4],  a  few  control  points  were  used  in 
order  to  ensure  an  underdetermined  system  of  equations 
and  obtain  an  exact  beam  pattern  at  those  few  points;  in 
our  technique,  the  entire  goal  psf  is  used  to  obtain  the  least 
squares  solution  of  an  overdetermined  system  of  equations. 
This  method  enables  excellent  control  of  the  system  psf, 
and  has  a  significant  impact  on  aperture  design  for  several 
applications  such  as  improved  DOF.  Simulation  results  for 
these  and  other  examples  are  described  in  an  accompany¬ 
ing  paper  [2]. 


V.  Conclusions 

The  MSSE  technique  is  a  general  beamforming  method 
that  can  be  used  to  design  apertures  for  specific  appli¬ 
cations.  It  enables  the  design  of  arbitrary  beam  profiles 
by  calculating  the  appropriate  optimum  aperture  weight¬ 
ings.  The  system  performance  is  optimized  because  the 
calculated  weightings  minimize  the  SSE  between  the  de¬ 
sired  and  achieved  system  responses.  The  algorithm  can 
be  implemented  readily  in  both  CW  and  broadband  sys¬ 
tems.  In  CW  systems,  the  receive  weights  can  be  imple¬ 
mented  through  apodization  and  time  delays  or  complex 
weights.  In  broadband  systems,  implementation  is  analo¬ 
gous  to  having  a  dynamic  FIR  filter  on  each  channel. 
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Abstract — In  the  first  part  of  this  work,  we  introduced 
the  minimum  sum  squared  error  (MSSE)  technique  of  ul¬ 
trasound  beamformer  design.  This  technique  enables  the 
optimal  design  of  apertures  to  achieve  arbitrary  system 
responses.  In  the  MSSE  technique,  aperture  weights  are 
calculated  and  applied  to  minimize  the  sum  squared  error 
(SSE)  between  the  desired  and  actual  system  responses. 
In  this  paper,  we  present  the  results  of  simulations  per¬ 
formed  to  illustrate  the  implementation  and  validity  of  the 
MSSE  technique.  Continuous  wave  (CW)  and  broadband 
simulations  are  presented  to  demonstrate  the  application 
of  the  MSSE  method  to  obtain  arbitrary  system  responses 
(within  fundamental  physical  limitations  of  the  system).  We 
also  describe  CW  and  broadband  simulations  that  imple¬ 
ment  the  MSSE  method  for  improved  conventional  depth 
of  field  (DOF)  and  for  improved  correlation  DOF  in  trans¬ 
lated  aperture  geometries.  Using  the  MSSE  technique,  we 
improved  the  conventional  DOF  by  more  than  200%  in  CW 
simulations  and  more  than  100%  in  broadband  simulations. 
The  correlation  DOF  in  translated  aperture  geometries  was 
improved  by  more  than  700%  in  both  CW  and  broadband 
simulations. 


I.  Introduction 

IN  an  accompanying  paper  [1],  we  described  the  mini¬ 
mum  sum  squared  error  (MSSE)  technique  of  aperture 
design.  This  technique  enables  the  design  of  ultrasound 
systems  for  arbitrary  system  responses  and  beam  profiles. 
This  is  done  by  calculating  and  applying  the  aperture 
weightings  that  optimally  produce  these  profiles  within 
fundamental  limits  imposed  by  wavelength.  In  the  MSSE 
design  technique,  the  system  point  spread  function  (psf)  is 
expressed  in  a  linear  algebra  formulation  as  a  function  of 
aperture  weights  and  a  propagation  function.  This  prop¬ 
agation  function  may  be  derived  theoretically  or  deter¬ 
mined  by  simulations  or  experiments.  The  sum  squared 
error  (SSE)  between  the  desired  (goal)  psf  and  the  real¬ 
ized  psf  is  then  minimized.  A  brief  review  of  the  major 
results  derived  in  [1]  for  both  continuous  wave  (CW)  and 
broadband  systems  is  provided  below. 
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A.  One-Way  Continuous  Wave  Formulation 

We  can  express  the  one-way  transmit  psf  at  some  range 
z  as  a  function  of  a  propagation  matrix  Sz  and  the  trans¬ 
mit  aperture  weightings  T.  Sz  may  include  factors  such 
as  limited  element  angular  response  [2]  and  frequency  de¬ 
pendent  attenuation.  This  relationship  is  expressed  in  a 
linear  algebra  formulation  as  a  matrix  multiplication.  Us- 
ing  [3],  we  can  .then  derive  the  least  squares  solution  for 
the  transmit  weightings  that  minimize  the  sum  squared  er¬ 
ror  (SSE)  between  the  system  psf  and  the  desired  or  goal 
psf.  Because  the  SSE  is  minimized,  application  of  these 
transmit  weightings  results  in  the  generation  of  a  one-way 
psf  that  is  optimally  similar  to  the  goal  psf.  The  transmit 
weightings  that  minimize  the  SSE  between  the  system  psf 
and  goal  psf  are  given  by: 

T  =  (Sf  Sz)~l  S? Ptz  =  Sf  Ptz,  (1) 

where  the  superscripts  H>  —  1,  and  #  denote  the  conjugate 
transpose,  matrix  inverse  and  pseudoinverse  operations, 
respectively,  and  Ptz  represents  the  goal  psf.  Sf  is  the 
pseudoinverse  of  Sz. 

B.  Two-Way  Continuous  Wave  Formulation 

We  also  can  derive  the  receive  aperture  weights  that 
would  minimize  the  SSE  between  the  two-way  system  psf 
and  a  two-way  desired  or  goal  psf.  The  two-way  psf  can 
be  expressed  using  the  well-known  radar  equation  [4].  The 
radar  equation  states  that  the  two-way  psf  is  the  product 
of  the  one-way  transmit  and  receive  psfe.  The  receive  psf 
is  expressed  as  a  function  of  a  set  of  receive  weights  and 
a  receive  propagation  matrix  in  a  similar  fashion  to  the 
transmit  psf  described  in  the  previous  section.  The  trans¬ 
mit  psf  is  diagonalized  to  eliminate  the  point-by-point  mul¬ 
tiplication  operation,  and  then  combined  with  the  receive 
propagation  matrix.  This  process  yields  an  expression  that 
formulates  the  two-way  psf  in  a  linear  algebra  formulation 
as  a  function  of  the  receive  weights  and  the  secondary  func¬ 
tion  derived  from  the  one-way  transmit  psf  and  the  receive 
propagation  matrix.  The  SSE  between  the  two-way  psf  and 
the  goal  psf  is  minimized  by  the  following  receive  weights: 

R  ~  {l*TzdS  ErzdtS )  I*TzdS  EtRz  —  &TzdS  ^THzi 

(2) 

where  PrxdS  is  the  secondary  function  described  above 
and  is  constructed  by  diagonalizing  Ptz  with  its  elements 


0885-30 10/$  10. 00  ©  2003  IEEE 


26 


IEEE  TRANSACTIONS  ON  ULTRASONICS,  FERROELECTRJCS,  AND  FREQUENCY  CONTROL,  VOL.  50,  NO.  1,  JANUARY  2003 


along  the  0th  diagonal,  and  then  multiplying  the  result¬ 
ing  matrix  with  the  propagation  function  Sz.  P^zdS  is  the 
pseudoinverse  of  PrzdS- 

C.  One -  Way  Broadband  Formulation 

The  above  CW  formulations  are  of  limited  use  in  broad¬ 
band  ultrasound  systems  that  utilize  finite  bandwidth 
pulses.  Because  CW  ultrasound  is  used  only  for  limited 
applications  such  as  CW  Doppler,  we  also  derived  one-way 
and  two-way  formulations  that  adapt  the  MSSE  technique 
for  broadband  systems.  The  one-way  psf  can  be  expressed 
as  the  product  of  a  propagation  matrix  Az  with  the  trans¬ 
mit  aperture  weights  T.  The  psf  is  a  function  of  lateral  po¬ 
sition  and  time,  and  the  aperture  weights  arc  a  function  of 
the  element  number  and  time.  The  propagation  matrix  de¬ 
pends  upon  the  excitation  pulse  and  the  impulse  responses 
of  the  elements  that  comprise  the  transmit  aperture.  Using 
[3],  the  transmit  weights  that  minimize  the  SSE  between 
the  system  psf  and  the  goal  psf  can  be  solved  for  and  are 
given  by: 

T  =  {ATZ  Az)~l  Atz  PTz  =  A*  PTz ,  (3) 

where  Prz  is  the  goal  or  desired  psf,  and  Af  is  the  pseu¬ 
doinverse  of  Az.  Implementation  of  these  weights  is  analo¬ 
gous  to  applying  a  dynamic  Finite  Impulse  Response  (FIR) 
filter  to  every  channel  and  summing  their  outputs. 

D.  Two-Way  Broadband  Formulation 

.  The  two-way  pulse-echo  psf  can  be  expressed  as  the 
product  of  a  two-way  propagation  matrix  Azz  with  the  re¬ 
ceive  aperture  weights  R.  Azz  depends  upon  the  excitation 
pulse,  the  transmit  and  receive  aperture  element  impulse 
responses,  and  the  transmit  aperture  weightings.  The  re¬ 
ceive  weightings  that  optimize  the  system  psf  are  given  by: 

R  =  (y£  ATZZ  PTRz  =  A*  PTRz,  (4) 

where  Ptrz  Is  the  goal  two-way  psf  and  Afz  is  the  pseu¬ 
doinverse  of  Azz . 

E.  Reduced  Computational  Cost 
Through  Symmetry  Relations 

The  design  of  weights  via  the  MSSE  technique  is  com¬ 
putationally  expensive  because  of  the  need  to  compute  the 
pseudoinverse  of  the  large  propagation  matrix  during  the 
calculation  of  the  aperture  weights.  In  order  to  reduce  com¬ 
putational  complexity,  we  exploit  the  lateral  symmetry  of 
the  transducer  apertures  and  the  psfe  [1],  Due  to  this  sym¬ 
metry,  it  is  sufficient  to  compute  only  half  of  the  weights, 
using  just  half  of  the  psf.  This  results  in  the  reduction 
of  the  computational  complexity  and  enables  more  effi¬ 
cient  computation  of  the  aperture  weights.  Note  that  the 
application  of  the  weights  has  a  relatively  low  computa¬ 
tional  cost  when  compared  to  the  design  of  the  weights. 


TABLE  I 

Parameters  Used  in  Continuous  Wave  (CW)  Simulations, 
Unless  Otherwise  Mentioned. 


Parameter 

Value 

Number  of  elements 

32 

Element  pitch 

135  pm 

Focus 

1.3  cm 

Lateral  window  over  which  the  psf  was  calculated 

180° 

psf  window  sampling  interval 

0.01° 

Ultrasonic  wave  propagation  speed 

1540  m/s 

Frequency 

10  MHz 

TABLE  II 

Parameters  Used  in  Broadband  Simulations,  Unless 
Otherwise  Mentioned. 

Parameter 

Value 

Number  of  elements 

32 

Element  pitch 

135  pm 

Focus 

1.3  cm 

Lateral  window  over  which  the  psf  was  calculated 

70  mm 

psf  window  sampling  interval 

50  pm 

Ultrasonic  wave  propagation  speed 

1540  ra/s 

Center  frequency 

10  MHz 

-6  dB  Bandwidth 

75% 

Ttemporal  sampling  of  pef 

120  MHz 

Temporal  spacing  of  weights  for  each  element 

25  ns 

The  weights  can  either  be  computed  and  stored  in  a  look¬ 
up  table  system  or  calculated  as  required.  They  can  then 
be  applied  as  conventional  FIR  filters  in  the  one-way  de¬ 
sign  case,  or  as  dynamic  shift  variant  filters  in  the  two-way 
design  case. 

All  the  formulas  presented  above  are  comprehensively 
derived  in  [l].  This  paper  describes  the  results  of  simula¬ 
tions  that  were  implemented  to  demonstrate  the  validity 
and  flexibility  of  the  MSSE  design  technique. 


II.  Simulations 

We  performed  two  sets  of  simulations  to  investigate  the 
performance  of  the  MSSE  technique.  The  first  set  was  in¬ 
tended  to  illustrate  the  implementation  of  the  MSSE  tech¬ 
nique  to  obtain  a  predetermined  system  psf  in  CW  and 
broadband  systems.  The  second  set  was  designed  to  imple¬ 
ment  the  technique  in  some  of  the  examples  of  application 
described  in  [1],  The  default  CW  system  parameters  are 
described  in  Table  I,  and  the  default  broadband  system 
parameters  are  described  in  Table  II.  Unless  stated  oth¬ 
erwise,  these  parameters  were  used  in  all  simulations.  All 
simulations  were  performed  in  Matlab  (The  Mathworks, 
Inc.,  Natick,  MA).  We  used  Field  II,  an  ultrasound  simula¬ 
tion  package  developed  by  Jensen  [5] ,  in  all  the  broadband 
simulations. 
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A.  Single  Range  Design  Examples 

1.  GW  One-Way  Design  Example  (Transmit  only):  For 
all  CW  simulations,  the  propagation  function  was  derived 
from  the  Rayleigh-Sommerfeld  equation  ([6,  pp.  46-50]). 
The  elements  were  treated  as  point  sources.  The  propaga¬ 
tion  function  for  the  field  at  an  arbitrary  point  p  due  to 
element  e  is: 

_  exp  (jkrep)  .  . 

Tep 

where  rcp  is  the  length  of  the  vector  pointing  from  element 
e  to  field  point  p,  and  k  =  ^  is  the  wave  number,  with  A 
being  the  acoustic  wavelength. 

We  used  the  parameters  described  in  Table  I  in  the  sim¬ 
ulations.  We  chose  the  goal  psf  to  be  a  Hann  window  [7,  pp. 
170-172]  of  width  32°  and  calculated  the  optimum  trans¬ 
mit  weights  using  (1).  No  transmit  apodization  other  than 
the  calculated  weights  was  used.  The  window  of  analysis 
was  ±90°.  Due  to  the  large  angles  of  interrogation,  we  in¬ 
cluded  the  obliquity  factor  in  these  simulations.  The  obliq¬ 
uity  factor  accounts  for  limited  element  angular  response 
and  was  included  in  the  propagation  function  as  follows: 

(t) 

Tep 

where  cos  (Oep)  is  the  obliquity  factor,  and  0ep  is  the  angle 
between  the  vector  pointing  from  element  e  to  field  point 
p ,  and  the  vector  normal  to  the  plane  of  the  transducer. 
Because  elements  were  modeled  as  point  sources,  the  sine 
directivity  pattern  expected  for  elements  of  finite  extent 
was  not  used. 

We  calculated  the  system  psf  after  applying  the  MSSE 
designed  transmit  weightings.  Figs.  1(a)  and  1(b)  show 
the  goal  psf  and  the  designed  system  psf,  respectively. 
Figs.  1(c)  and  1(d)  display  the  magnitude  and  phase  of 
the  calculated  transmit  aperture  weights. 

Performance  of  the  MSSE  technique.  In  an  effort  to 
explore  the  performance  of  the  algorithm,  we  made  the 
Hann  window  goal  psf  increasingly  narrower.  We  progres¬ 
sively  reduced  the  goal  psf  width  to  16°,  8°,  and  4°.  The 
Rayleigh  resolution  limit  for  the  geometry  that  was  used 
is  around  2.04°.  Figs.  1(e)  and  1(f)  show  the  goal  and  gen¬ 
erated  system  psfe,  respectively,  when  the  goal  psf  was  16° 
wide.  Figs.  1(g)  and  1(h)  show  the  magnitude  and  phase 
of  the  calculated  weights,  respectively,  for  the  16°  goal  psf. 
Figs.  l(i)-l(l)  display  the  same  information  for  a  goal  psf  of 
width  8°,  and  Figs.  l(m)~l(p)  depict  the  results  obtained 
when  the  goal  psf  had  a  width  of  4°. 

Effects  of  the  size  of  the  window  of  analysis .  The  MSSE 
algorithm  optimizes  the  system  psf  only  within  the  win¬ 
dow  of  analysis.  Effects  that  occur  outside  this  window  are 
ignored,  thereby  introducing  potential  artifacts  in  the  ul¬ 
trasonic  field  generated  outside  the  analysis  window.  We 
performed  simulations  to  investigate  the  effects  of  varying 
this  window’s  size  on  the  ultrasonic  field  outside  the  win¬ 
dow.  The  goal  psf  was  a  Hann  window  that  was  6°  wide. 


The  MSSE  technique  was  implemented  for  analysis  win¬ 
dow  sizes  of  ±15°,  ±30°,  ±45°,  ±60°,  ±75°,  and  ±90°, 
sampled  every  0.01°.  The  system  psf  was  then  computed 
over  a  ±90°  window  using  the  calculated  transmit  weights. 

Fig.  2(a)  depicts  the  ideal  psf  used  for  a  ±15a  window  of 
analysis.  Fig.  2(b)  illustrates  the  psf  obtained  after  imple¬ 
mentation  of  the  MSSE  algorithm,  and  Fig.  2(c)  shows  the 
obtained  psf  on  a  logarithmic  scale  after  being  normalized 
to  the  peak  mainlobe  level.  Figs.  2(d)-2(f),  2(g)-2(i),  2(j)- 
2(1),  2(m)-2(o),  and  2(p)-2(r)  depict  the  same  information 
for  window  analysis  sizes  of  ±30°,  ±45°,  ±60°,  ±75°,  and 
±90°,  respectively.  The  extent  of  the  window  of  analysis 
is  shown  by  dotted  lines  in  each  plot. 

Effect  of  errors  in  the  assumed  ultrasonic  wave  propa¬ 
gation  speed .  The  adverse  effects  of  errors  in  the  assumed 
wave  propagation  speed  on  the  response  of  an  ultrasound 
system  are  well-known  [8],  Because  the  MSSE  technique 
uses  dynamic  shift-variant  aperture  weights,  errors  in  the 
assumed  wave  speed  are  an  important  concern.  There¬ 
fore,  we  implemented  simulations  in  which  the  actual  wave 
speed  was  underestimated  by  25  m/s,  50  m/s,  and  75  m/s, 
and  then  overestimated  by  25  m/s,  50  m/s,  and  75  m/s. 
The  goal  psf  was  a  6°  wide  Hann  window.  The  window 
of  analysis  was  ±90°,  and  it  was  sampled  every  0.01°.  We 
compared  the  psfe  obtained  in  these  simulations  to  the  psf 
obtained  when  the  assumed  speed  was  correct,  which  is 
shown  in  Fig.  3(a).  Figs.  3(c),  3(g),  and  3(k)  depict  the 
psfe  obtained  when  the  assumed  sound  speed  was  under¬ 
estimated  by  25  m/s,  50  m/s,  and  75  m/s,  respectively. 
Figs.  3(e),  3(i),  and  3(m)  display  the  psfe  obtained  when 
the  assumed  speed  was  overestimated  by  25  m/s,  50  m/s, 
and  75  m/s,  respectively.  Fig.  3(b)  depicts  the  obtained 
psf  error  magnitude  when  the  assumed  speed  was  correct. 
Figs.  3(d),  3(h),  and  3(1)  show  the  psf  error  magnitude 
when  the  speed  was  underestimated,  and  Figs.  3(f),  3(j), 
and  3(n)  display  the  psf  error  magnitude  when  the  speed 
was  overestimated.  The  presented  psfe  were  computed  at 
the  intended  transmit  foci,  not  the  shifted  foci.  The  ob¬ 
jective  of  these  simulations  was  to  investigate  the  degra¬ 
dation  in  the  performance  of  the  MSSE  algorithm  that  is 
produced  by  an  incorrect,  assumed  propagation  speed  at  a 
particular  location  of  interest.  Simulations  were  also  per¬ 
formed  at  the  shifted  focus.  Although  the  results  are  not 
presented  here,  they  were  qualitatively  similar. 

2.  CW  Two-Way  Design  Example  (Transmit-Receive): 

A  Hann  window  of  width  4°  was  chosen  to  be  the  goal 
two-way  psf.  We  implemented  the  MSSE  technique  by  cal¬ 
culating  and  applying  optimum  receive  aperture  weights 
using  (2).  No  apodization  was  used  on  the  transmit  aper¬ 
ture.  Table  I  lists  the  parameters  used  in  the  simulation. 
We  then  computed  the  two-way  system  psf  using  the  calcu¬ 
lated  receive  weights.  Figs.  4(a)  and  4(b)  depict  the  goal 
psf  and  the  achieved  system  psf.  Fig.  4(c)  displays  the 
magnitude  of  the  error  between  the  goal  and  system  psfe. 
Figs.  4(d)  and  4(e)  display  the  magnitude  and  phase  of  the 
calculated  receive  aperture  weights,  respectively. 
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Goalpsfs  Obtained  psfs  Calculated  Tx  weights  Calculated  Tx  weights 

(Magnitude)  (Phase) 


00  00  (c)  (d) 


©  0)  (k)  (1) 


(m)  (n)  (o)  (p) 


Fig.  1.  One-way  CW  design  example  for  goal  psfs  of  width  32°,  16°,  8°,  and  4°.  The  first  column  depicts  the  goal  psfs;  the  second  column 
shows  the  obtained  psfs,  both  as  a  function  of  lateral  position.  Results  for  goal  psf  widths  of  32°,  16°,  8°,  and  4°,  respectively,  are  shown 
from  top  to  bottom.  The  third  and  fourth  columns  illustrate  the  magnitude  and  phase  of  the  calculated  transmit  weights,  respectively  as  a 
function  of  element  number.  It  can  be  seen  that  the  obtained  psfs  are  qualitatively  similar  to  the  goal  psfs  in  all  the  simulations.  The  profile 
of  the  magnitude  of  the  calculated  weights  shown  in  (c),  (g),  (k),  and  (o)  expectedly  becomes  wider  as  the  goal  psf  becomes  narrower 
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Fig.  2.  Effects  of  the  size  of  the  window  of  analysis.  The  first  column  depicts  the  goal  psfs  while  the  second  column  shows  the  obtained 
psfs,  respectively,  both  as  a  function  of  lateral  position,  for  different  analysis  window  sizes.  Results  for  analysis  window  sizes  of  ±15°,  ±30° , 
±45°,  ±60°,  ±75°,  and  ±90°are  shown  from  top  to  bottom.  The  third  column  shows  the  obtained  psfs  on  a  logarithmic  scale  after  being 
normalized  to  the  peak  mainlobe  level.  The  edges  of  the  window  of  analysis  are  indicated  by  the  dotted  lines  in  each  plot.  It  can  be  seen 
dearly  that  the  size  of  the  analysis  window  plays  a  critical  role  in  the  performance  of  the  MSSE  algorithm.  As  can  be  seen  from  plots  (b), 
(c),  (e),  and  (f),  the  obtained  psf  has  large  grating  lobes  outside  the  analysis  window  if  it  is  too  small.  An  analysis  window  of  at  least  ±45° 
is  required. 


II] 


IEEE  TRANSACTIONS  ON  ULTRASONICS,  FERROELECTRICS,  AND  FREQUENCY  CONTROL,  VOL.  50,  NO.  1,  JANUARY  2003 


Obtained  psf 


Error 


Lateral  position  (degrees)  Lateral  position  (degrees) 


(a) 


(b) 


Obtained  psfs 


Errors 


| 

'c 

CD 

«c 


Il  L —  crisis  11 

i 

5 

L 

Obtained  psfs 

f—  C=1565  II 


Errors 


(c) 


(i) 


0) 


IS 


r  j  —  G-146S  II 

Lateral  position  (degrees) 

(k) 


-60  -30  0  30  60  90 

Lateral  position  (degrees) 


*60  -30  0  30  60  90 

Lateral  position  (degrees) 


(i) 


(m) 


Lateral  position  (degrees) 

(n) 


Fig.  3.  Effect  of  errors  in  the  assumed  wave  propagation  speed.  The  plots  show  the  obtained  psfs  and  psf  errors  when  the  assumed  wave 
speed  is  correct,  underestimated  by  25  m/s,  50  m/s,  and  75  m/s,  and  overestimated  by  25  m/s,  50  m/s,  and  75  m/s.  (a)  depicts  the  obtained 
psf  with  the  correct  wave  speed,  and  (b)  illustrates  the  obtained  psf  error.  Below,  the  first  and  second  columns  depict  the  obtained  psfs  and 
psf  error  magnitudes  as  a  function  of  angle  when  the  speed  was  underestimated  by  25  ra/s,  50  m/s,  and  75  m/s,  respectively.  The  third 
and  fourth  columns  depict  the  obtained  psfs  and  psf  error  magnitudes  as  a  function  of  angle  when  the  speed  was  overestimated  by  25  m/s, 
50  m/s,  and  75  m/s  respectively.  The  performance  of  the  MSSE  technique  was  slightly  degraded  when  the  assumed  wave  propagation  speed 
was  incorrect,  but  the  overall  errors  were  reasonable,  showing  that  the  MSSE  technique  is  stable. 
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Fig.  4.  Two-way  CW  design  example,  (a)  depicts  the  goal  psf  as  a 
function  of  lateral  position,  and  (b)  depicts  the  obtained  psf  as  a 
function  of  lateral  position,  (c)  illustrates  the  error  between  the  goal 
and  the  obtained  psfs  as  a  function  of  lateral  position,  (d)  shows  the 
magnitude  and  (e)  shows  the  phase  of  the  calculated  receive  weights, 
respectively,  as  a  function  of  element  number. 

5.  Broadband  One-Way  Example  (Transmit  only):  We 
used  the  ultrasound  simulation  package  Field  II,  developed 
by  Jensen  [5] ,  for  all  broadband  simulations.  All  broadband 
simulations  took  advantage  of  system  symmetry  to  reduce 
computational  complexity  [1],  as  previously  described.  The 
propagation  matrix,  Az,  was  a  four-dimensional  function. 
For  an  aperture  of  N  elements,  every  term  of  Az  was  of  the 
form  if  determined  the  field  due  to  the  weight¬ 

ing  at  time  sample  j  applied  to  elements  i  and  N  —  i,  at 
time  sample  l  at  lateral  field  point  k.  It  was  constructed  us¬ 
ing  dual  element  spatial  impulse  responses  that  described 
the  contributions  of  a  selected  pair  of  elements  at  each  field 
point  at  each  sampled  time  point.  The  dual  element  re¬ 
sponses  were  determined  by  transmitting  only  on  selected 
pairs  of  elements.  These  were  then  used  to  construct  the 


propagation  function.  The  goal  psf  was  generated  by  axi¬ 
ally  weighting  a  sinusoidal  signal  by  a  Hann  window,  and 
multiplying  the  result  by  a  lateral  Hann  window.  Note 
that  this  goal  psf  is  quite  challenging  because  it  lacks  the 
wavefront  curvature  that  would  normally  be  expected. 

In  all  broadband  simulations,  we  downsampled  the 
propagation  matrix  by  a  factor  of  3.  This  had  the  effect  of 
reducing  the  upper  cut-off  frequency  in  the  frequency  re¬ 
sponse  of  the  FIR  filter  formed  by  the  calculated  weights. 
The  temporal  sampling  rate  of  the  psf  was  120  MHz.  The 
upper  cut-off  frequency,  therefore,  was  reduced  from  half 
the  temporal  sampling  rate  of  the  psf  (60  MHz)  to  20  MHz. 
This  rate  still  provided  adequate  sampling  given  an  input 
pulse  with  a  center  frequency  of  10  MHz  and  a  —6  dB  rel¬ 
ative  bandwidth  of  75%.  Table  II  lists  the  parameters  used 
in  the  simulation. 

We  then  calculated  the  optimum  transmit  weights  us¬ 
ing  the  goal  psf  and  the  downsampled  propagation  matrix. 
Figs.  5(a)  and  5{b)  show  the  goal  and  the  achieved  system 
psfs,  respectively,  both  as  a  function  of  lateral  position  and 
time.  We  also  envelope-detected  the  psfs  using  the  Hilbert 
transform  ([9,  pp.  359-367]),  and  then  peak-detected  them 
in  the  time  dimension  to  generate  beam  profiles.  Figs.  5(c) 
and  5(d)  show  these  goal  and  achieved  system  psf  profiles 
as  a  function  of  lateral  position.  Figs.  5(e)  and  5(f)  dis¬ 
play  the  calculated  transmit  weights  as  a  function  of  the 
element  number  and  time,  and  the  magnitude  of  the  error 
as  a  function  of  lateral  position  and  time,  respectively.  The 
error  was  calculated  by  computing  the  difference  between 
the  goal  and  system  psfs.  Note  that  the  weights  do  not 
include  geometric  focal  delays.  The  delays  were  applied 
separately  prior  to  the  application  of  the  weights.  The  re¬ 
sult  obtained  by  convolution  of  the  calculated  weights  with 
the  transmit  pulse  is  presented  in  Fig.  5(g). 

4.  Broadband  Two-Way  Example  (Transmit- Receive): 
The  goal  psf  was  constructed  in  a  similar  manner  to 
the  one-way  example,  except  for  a  scaling  factor  that  ac¬ 
counted  for  the  reduction  in  magnitude  due  to  two-way 
propagation.  The  propagation  matrix  Azz,  however,  was 
different  from  the  one-way  case.  We  constructed  the  prop¬ 
agation  functions  by  using  the  entire  transmit  aperture 
and  receiving  only  on  selected  pairs  of  elements.  Transmit 
focal  delays  were  applied,  but  no  transmit  apodization  was 
used.  All  the  other  parameters  are  listed  in  Table  II.  We 
then  downsampled  the  resulting  propagation  matrix  by  a 
factor  of  3  and  used  it  to  calculate  the  receive  weights 
that  minimize  the  SSE.  Figs.  6(a)  and  6(b)  display  the 
goal  and  achieved  system  psfs,  respectively.  As  in  the  one¬ 
way  case,  we  envelope-detected  and  peak-detected  the  psfs. 
The  goal  and  achieved  system  psf  profiles  are  displayed  in 
Figs.  6(c)  and  6(d),  respectively.  Figs.  6(e)  and  6(f)  show 
the  calculated  receive  weights  and  the  psf  error  magnitude 
respectively. 

B.  Applications 

1.  Enhanced  Depth  of  Field:  The  DOF  of  an  ultrasound 
system  is  the  axial  region  within  which  the  system  is  said 
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Pig.  5.  One-way  broadband  design  example,  (a)  depicts  the  goal  psf  as  a  function  of  lateral  position  and  time,  and  (b)  depicts  the  obtained 
psf  as  a  function  of  lateral  position  and  time,  (c)  shows  the  goal  psf  beam  profile  as  a  function  of  lateral  position,  and  (d)  shows  the  obtained 
psf  beam  profile  as  a  function  of  lateral  position.  These  profiles  were  obtained  by  envelope  detection  using  the  Hilbert  transform  (9,  pp. 
359-367],  followed  by  peak  detection  of  the  psf  along  the  time  dimension,  (e)  shows  the  calculated  receive  weights  as  a  function  of  element 
number  and  time,  and  (f)  illustrates  the  error  between  the  goal  and  the  obtained  psfe  as  a  function  of  lateral  position  and  time,  (g)  depicts 
the  results  obtained  by  convolving  the  calculated  transmit  weights  with  the  transmit  pulse.  It  can  be  seen  from  (a)  and  (b)  that  the  goal 
and  obtained  psfs  are  very  similar.  The  beam  profile  of  the  obtained  psf  in  (d)  shows  high  frequency  errors,  which  are  also  noticeable  in 
the  error  image  in  (f).  These  errors  arc  caused  by  the  inadequately  low  temporal  sampling  rate  used  in  our  simulations  in  Field  II  due  to 
computational  limitations,  and  should  be  eliminated  with  a  higher  sampling  rate. 
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Fig.  6.  Two-way  broadband  design  example,  (a)  depicts  the  goal  psf  as  a  function  of  lateral  position  and  time,  and  (b)  depicts  the  obtained 
psf  as  a  function  of  lateral  position  and  time,  (c)  shows  the  goal  psf  beam  profile  as  a  function  of  lateral  position,  and  (d)  shows  the  obtained 
psf  beam  profile  as  a  function  of  lateral  position.  These  profiles  were  obtained  by  envelope  detection  using  the  Hilbert  transform  (9,  pp. 
359-367],  followed  by  peak  detection  of  the  psf  along  the  time  dimension,  (e)  shows  the  calculated  receive  weights  as  a  function  of  element 
number  and  time,  and  (f)  illustrates  the  error  between  the  goal  and  the  obtained  psfs  as  a  function  of  lateral  position  and  time.  It  can  be 
seen  from  (a)  and  (b)  that  the  goal  and  obtained  psfs  are  very  similar.  As  in  the  one-way  simulations,  the  beam  profile  of  the  obtained 
psf  in  (d)  shows  high  frequency  errors,  which  are  also  noticeable  in  the  error  image  in  (f).  These  errors  are  caused  by  the  inadequately 
low  temporal  sampling  rate  used  in  our  simulations  in  Field  II  due  to  computational  limitations,  and  should  be  eliminated  with  a  higher 
sampling  rate. 
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to  be  in  focus,  or  the  axial  region  within  which  the  sys¬ 
tem  psf  satisfies  a  chosen  criterion.  The  system  psf  should 
ideally  remain  similar  to  the  psf  at  the  focus  for  as  large 
an  axial  range  as  possible.  Application  of  the  MSSE  tech¬ 
nique  using  the  psf  at  the  focus  as  the  goal  psf  at  each 
range  should  theoretically  yield  weights  that  minimize  the 
SSE  between  the  goal  psf  and  the  psf  at  the  range  of  in¬ 
terest.  This,  in  turn,  should  maximize  the  DOF.  With  this 
rationale,  we  applied  the  two-way  MSSE  technique  in  both 
CW  and  broadband  simulations  in  an  attempt  to  improve 
the  DOF  of  the  system.  The  technique  was  applied  at  ev¬ 
ery  sampled  axial  range  point.  The  goal  or  desired  psf  at 
each  range  point  was  the  psf  computed  at  the  focus. 

In  both  CW  and  broadband  simulations,  we  generated 
the  goal  psf  using  10-element  transmit  and  receive  aper¬ 
tures,  while  we  used  a  16-element  transmit  and  32-element 
receive  aperture  for  the  actual  system  psf.  i.c.  weights 
were  calculated  for  32  receive  elements.  The  system  was 
focused  at  6.5  mm.  A  Hann  window  was  used  as  a  trans¬ 
mit  apodization  function  in  all  simulations,  and  as  receive 
apodization  during  the  generation  of  the  goal  psf.  We  per¬ 
formed  CW  simulations  to  maximize  the  DOF  over  an 
axial  window  from  0.1  mm  to  50  mm  that  was  sampled 
every  0.1  mm.  The  lateral  window  over  which  the  CW  psf 
was  calculated  was  sampled  at  5  fim  to  enable  accurate 
estimation  of  the  full  width  at  half  maximum  (FWHM) 
of  the  mainlobe.  We  also  implemented  broadband  simula¬ 
tions  over  an  axial  window  from  0.5  mm  to  32.5  mm  that 
was  sampled  every  2  mm.  All  other  parameters  are  listed 
in  Tables  I  and  II. 

Current  ultrasound  systems  attempt  to  improve  the 
DOF  by  using  dynamic  apodization  and  dynamic  receive 
focusing  [10],  and  we  used  these  in  control  simulations  to 
establish  a  basis  for  comparison  with  the  results  obtained 
using  the  MSSE  technique.  We  allowed  the  apodization 
profile  to  grow  as  a  function  of  range  assuming  an  infinitely 
large  aperture,  and  used  the  central  portion  of  the  profile 
to  implement  dynamic  apodization  on  the  32  receive  ele¬ 
ments. 

In  broadband  simulations,  the  temporal  sampling  was 
120  MHz,  and  the  propagation  matrix  was  downsampled 
by  a  factor  of  3  to  limit  the  frequency  response  of  the  cal¬ 
culated  weights.  In  order  to  perform  accurate  analysis,  we 
interpolated  the  broadband  psfs  by  a  factor  of  10  in  both 
the  temporal  and  lateral  spatial  dimensions  using  cubic 
spline  interpolation. 

We  calculated  correlation  coefficients  to  compare  the  psf 
at  each  range  of  interrogation  with  the  psf  at  the  focus. 
The  correlation  coefficients  were  computed  for  complex 
data  in  CW  simulations  and  for  real  data  in  broadband 
simulations.  We  also  calculated  the  FWHM  of  the  main- 
lobe  for  each  psf.  Figs.  7(a)  and  7(b)  depict  the  correla- 
tion  coefficients  and  the  mainlobe  FWTIM  obtained  in  CW 
simulations.  The  dotted  lines  indicate  the  transmit  focus. 
These  were  obtained  at  each  range  of  interrogation  for  the 
control  case  and  the  case  when  the  MSSE  technique  was 
applied.  Figs.  7(c)  and  7(d)  display  the  correlation  coeffi¬ 
cients  and  the  mainlobe  FWHM  for  broadband  simulations 


of  the  control  and  the  MSSE  technique  cases.  A  more  qual¬ 
itative  assessment  of  the  efficacy  of  the  MSSE  technique 
can  be  made  using  Figs.  7(c)  and  7(f).  Both  these  im¬ 
ages  were  constructed  by  superimposing  CW  psfs  at  mul¬ 
tiple  ranges  that  were  normalized  to  their  peak  mainlobe 
levels.  The  images  are  shown  on  a  logarithmic  scale  and 
are  a  function  of  range  and  lateral  position.  Fig.  7(e)  was 
constructed  using  psfs  obtained  in  the  control  simulation. 
Fig.  7(f)  was  generated  using  psfs  from  the  simulation  in¬ 
volving  the  MSSE  technique. 

2.  Increased  Correlation  Depth  of  Field  in  Translated 
Aperture  Geometries:  The  translating  apertures  algo¬ 
rithm  (TAA)  is  a  technique  that  is  used  in  phase  aber¬ 
ration  correction  [11],  [12]  and  to  obtain  accurate  angular 
scatter  data  [13].  Data  are  acquired  using  two  system  ge¬ 
ometries,  the  backscatter  geometry  in  which  the  transmit 
and  receive  apertures  are  coincident,  and  the  angular  scat¬ 
ter  geometry  in  which  the  apertures  are  translated  by  equal 
distances  in  opposite  directions.  Despite  its  many  advan¬ 
tages,  use  of  the  TAA  results  in  a  dramatically  reduced 
DOF.  This  is  due  to  the  separation  of  the  transmit  and 
receive  apertures,  which  causes  only  a  limited  interference 
of  the  transmit  and  receive  beams.  In  order  to  fully  utilize 
the  TAA,  a  high  correlation  between  the  psfs  of  the  trans¬ 
lated  and  the  nontranslated  geometries  is  required  over  an 
extended  axial  region. 

We  performed  CW  and  broadband  simulations  in  an  ef¬ 
fort  to  investigate  the  ability  of  the  MSSE  technique  to 
improve  the  correlation  DOF  in  the  TAA.  We  applied  the 
MSSE  technique  to  a  translated  apertures  system  with  the 
goal  psf  being  that  generated  by  the  same  system  without 
any  translation  of  the  apertures  at  the  same  range.  Both 
transmit  and  receive  apertures  comprised  16  elements.  In 
the  translated  geometry  case,  the  apertures  were  trans¬ 
lated  by  8  elements  in  opposite  directions.  The  system  was 
focused  at  6.5  mm.  As  in  the  enhanced  conventional  DOF 
case,  we  performed  CW  simulations  over  an  axial  window 
from  0.1  mra  to  50  mm  that  was  sampled  every  0.1  mm. 
We  also  implemented  broadband  simulations  over  an  axial 
window  from  0.5  mm  to  32.5  mra  that  was  sampled  ev¬ 
ery  2  mm.  A  Hann  window  was  used  as  apodization  on 
the  transmit  apertures  and  on  the  receive  aperture  dur¬ 
ing  the  generation  of  the  goal  psf.  The  temporal  sampling 
in  broadband  simulations  was  120  MHz.  We  downsampled 
the  broadband  propagation  matrix  by  a  factor  of  3,  and 
interpolated  the  system  psfs  by  a  factor  of  10  in  both  the 
time  and  lateral  space  dimensions  prior  to  analysis.  Other 
parameters  were  consistent  with  Tables  I  and  II. 

Correlation  coefficients  were  calculated  at  each  range  by 
correlating  the  psf  obtained  after  application  of  the  MSSE 
technique  in  the  translated  geometry  with  the  psf  obtained 
using  the  nontranslated  geometry.  As  a  comparison,  we 
also  performed  a  control  simulation  in  which  we  calcu¬ 
lated  correlation  coefficients  by  correlating  the  translated 
and  nontranslated  geometry  psfs,  without  the  application 
of  the  MSSE  technique.  Figs.  8(a)  and  8(b)  show  the  corre¬ 
lation  curves  obtained  in  CW  and  broadband  simulations, 
respectively.  The  dotted  lines  indicate  the  transmit  focus. 
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Fig.  7.  Application  of  the  MSSE  technique  for  enhanced  depth  of  field  (DOF),  (a)  shows  CW  correlation  coefficients,  calculated  by  correlating 
the  psf  at  each  range  of  interrogation  with  the  psf  at  the  focus,  as  a  function  of  range,  (b)  shows  the  FWHM  of  the  mainlobe  obtained  in 
CW  simulations  as  a  function  of  range,  (c)  shows  broadband  correlation  coefficients,  and  (d)  shows  the  FWHM  of  the  mainlobe  obtained  in 
broadband  simulations,  both  as  a  function  of  range.  Dynamic  apodization  and  dynamic  receive  focusing  were  used  in  the  control  simulations. 
The  dotted  lines  indicate  the  transmit  focus,  (e)  depicts  CW  psf  images  generated  using  control  psfs,  and  (f)  shows  psfs  obtained  by  applying 
the  MSSE  technique.  Both  images  were  formed  from  normalized  and  logarithmically  compressed  lateral  psfs  at  multiple  ranges.  Each  column 
in  the  images  consists  of  the  psf  at  a  single  range.  In  both  CW  and  broadband  simulations,  use  of  the  MSSE  technique  resulted  in  an  increased 
DOF.  This  is  apparent  from  the  correlation  coefficients  in  (a)  and  (c),  in  which  the  coefficients  obtained  from  the  MSSE  simulations  are  higher 
than  those  obtained  from  control  simulations.  It  also  can  be  seen  from  (b)  and  (c)  that  the  mainlobe  width  over  the  span  of  interrogated 
axial  ranges  is  maintained  much  better  in  MSSE  simulations  than  in  control  simulations.  A  qualitative  comparison  between  the  MSSE 
technique  and  control  simulations  can  be  made  from  (e)  and  (f). 


(<0  (d)  (e) 

Fig.  8.  Results  of  simulations  for  increased  correlation  depth  of  field  (DOF)  in  translated  aperture  geometries,  (a)  depicts  the  correlation 
coefficients  obtained  in  CW  simulations,  calculated  by  correlating  the  non- translated  geometry  psf  at  each  range  with  the  translated  geometry 
'  psf  at  the  same  range,  as  a  function  of  range,  (b)  shows  the  correlation  coefficients  obtained  in  broadband  simulations  as  a  function  of  range. 
The  dotted  lines  indicate  the  transmit  focus.  At  each  range,  the  MSSE  technique  was  applied  in  the  translated  geometry  case  with  the 
non-translated  geometry  psf  as  the  goal  psf.  In  control  simulations,  the  translated  and  non-translatcd  geometry  psfs  were  obtained  without 
application  of  the  MSSE  technique,  (c)  depicts  the  goal  psf  (which  is  the  psf  obtained  with  the  nontranslated  geometry),  the  translated 
geometry  control  psf,  and  the  translated  geometry  psf  obtained  using  the  MSSE  technique  at  an  axial  range  of  4  mm.  The  same  information 
at  the  transmit  focus  (6.5  mm)  and  at  9  mm  is  depicted  in  (d)  and  (e),  respectively.  It  is  apparent  that  the  correlation  coefficients  obtained 
from  the  MSSE  simulations  are  higher  than  those  obtained  from  control  simulations,  and  the  MSSE  technique  significantly  improves  the 
.  DOF.  It  can  be  seen  from  (c),  (d),  and  (e)  that  the  MSSE  algorithm  preferentially  optimizes  the  mainlobe  at  the  cost  of  increased  sidelobe 
levels.  If  it  is  important  to  maintain  low  sidelobe  levels,  the  weighted  MSSE  technique  described  in  [1]  can  be  used  to  achieve  a  compromise 
between  optimizing  the  mainlobe  and  controlling  sidelobe  levels. 

Fig.  8(c)  depicts  the  goal  psf  (the  psf  obtained  with  the  seen  from  Fig.  1(b)  that  the  designed  system  psf  closely 
nontranslated  geometry),  the  translated  geometry  control  approximates  it.  The  magnitude  and  phase  of  the  calcu- 
.  psf,  and  the  translated  geometry  psf  obtained  using  the  lated  transmit  weights  can  be  seen  in  Figs.  1(c)  and  1(d). 
MSSE  technique  at  an  axial  range  of  4  mm.  Figs.  8(d)  and  We  use  two  error  metrics  throughout  this  section.  We  use 

8(e)  depict  the  same  information  at  the  transmit  focus  the  ratio  of  the  root  mean  square  (rms)  error  between  the 

(6.5  mm),  and  at  9  mm.  psfs  to  the  peak  magnitude  of  the  goal  psf  as  a  metric  of 

the  error  in  our  desired  psf,  and  refer  to  it  as  the  relative 
rms  error.  We  also  use  the  ratio  of  the  peak  error  magni- 
III.  DISCUSSION  tude  to  the  peak  magnitude  of  the  goal  psf,  and  refer  to 

it  as  the  relative  peak  error.  The  relative  rms  error  in  the 

Figs.  l(a)-l(d)  demonstrate  the  use  of  the  MSSE  tech-  one-way  CW  simulation  was  0.029,  and  the  relative  peak 
nique  in  the  most  basic  system  configuration  that  was  sim-  error  was  0.066. 

ulated,  i.e.,  the  one-way  CW  system.  The  goal  psf  was  a  Figs.  l(a)-l(p)  illustrate  the  flexibility  of  the  MSSE 
32°-wide  Hann  window  as  shown  in  Fig.  1(a).  It  can  be  technique.  In  these  simulations,  the  goal  psf  width  was  pro- 
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gressively  reduced  to  16°,  8°,  and  4°  in  order  to  investigate 
the  performance  of  the  technique.  The  relative  rms  errors 
obtained  in  the  three  simulations  were  0.005,  0.0052,  and 
0.0295,  and  the  relative  peak  errors  were  0.012,  0.023,  and 
0.214,  respectively.  The  error  in  obtaining  the  16°  and  8° 
wide  goal  psfe  is  much  lower  than  the  error  for  the  3 2° and 
4°  wide  goal  psfs.  However,  it  can  be  seen  from  Fig.  1  that 
the  goal  and  obtained  psfs  are  qualitatively  quite  simi¬ 
lar.  An  examination  of  the  magnitude  of  the  calculated 
weights  in  Figs.  1(c),  1(g),  l(k),  and  l(o)  shows  that  the 
apodization  profile  expectedly  becomes  wider  as  the  goal 
psf  becomes  narrower. 

The  effects  of  varying  the  size  of  the  window  of  anal¬ 
ysis  are  shown  in  Fig.  2.  It  can  be  seen  from  Figs.  2(b) 
and  2(c)  that,  when  a  ±15°  analysis  window  was  used,  the 
ultrasonic  field  outside  the  window  was  unstable  and  unde¬ 
sirable,  exhibiting  erratic  behavior  with  very  large  grating 
lobes.  The  dotted  lines  indicate  the  extent  of  the  window 
of  analysis,  and  it  can  be  seen  that  the  grating  lobes  occur 
just  outside  the  window.  Results. improved  when  the  win¬ 
dow  size  was  progressively  increased  in  steps  of  ±15°  to 
±90°,  as  can  be  seen  from  Figs.  2{d)-2(r).  The  grating  lobe 
magnitudes  were  progressively  reduced  as  the  window  size 
was  increased.  Beyond  a  window  size  of  ±45°,  however,  the 
improvement  outside  the  analysis  window  was  negligible. 
In  the  ±15°  and  ±30°  window  cases,  it  can  be  seen  that 
the  grating  lobes  occur  immediately  outside  the  window. 
These  grating  lobes,  caused  by  inadequate  analysis  window 
size,  disappear  for  a  ±45°  window  of  analysis.  Therefore,  it 
can  be  seen  that  the  size  of  the  window  of  analysis  signif¬ 
icantly  impacts  the  obtained  ultrasonic  field,  and  it  must 
be  carefully  chosen  to  suit  the  application.  Ideally  it  would 
always  cover  ±90°,  although  computational  and  memory 
requirements  may  limit  the  practical  range. 

We  also  investigated  the  effect  of  errors  in  the  assumed 
speed  of  sound.  Simulation  results  are  displayed  in  Fig.  3. 
Fig.  3(b)  shows  the  psf  error  when  the  assumed  speed  of 
sound  was  correct.  In  this  case,  the  relative  rms  error  was 
0.012  and  the  relative  peak  error  was  0.054.  As  the  wave 
speed  was  underestimated  by  25  m/s,  50  m/s,  and  75  m/s, 
the  relative  rms  and  peak  errors  rose  to  0.022  and  0.062,  re¬ 
spectively,  This  is  shown  in  Figs.  3(d),  3(h)  and  3(1).  When 
the  wave  speed  was  overestimated  by  25m/s,  50  m/s,  and 
75  m/s,  the  relative  rms  and  peak  errors  rose  to  0.015 
and  0.102,  respectively,  as  can  be  seen  in  Figs.  3(f),  3(j) 
and  3(n).  The  observed  errors  were  reasonable,  despite  the 
extremely  wide  range  considered.  Therefore,  the  MSSE  al¬ 
gorithm  is  stable  in  the  sense  that  small  errors  in  the  as¬ 
sumed  wave  propagation  speed  do  not  appear  to  result  in 
a  significant  degradation  of  performance. 

Fig.  4  demonstrates  the  use  of  the  MSSE  technique  in 
the  design  of  two-way  system  responses.  It  can  be  seen 
from  Figs.  4(a)  and  4(b)  that  the  system  psf  obtained  by 
the  MSSE  technique  is  very  similar  to  the  goal  psf.  The  re¬ 
sulting  relative  rms  error  was  0.002,  and  the  relative  peak 
error  was  0.022.  Note  that  the  errors  obtained  in  the  one¬ 
way  and  two-way  simulations  cannot  be  compared  directly 
because  the  transmit  psf  was  predetermined  in  the  two¬ 


way  simulations,  and  the  two-way  psf  was  optimized  using 
only  the  receive  weights.  Therefore,  there  is  inherently  less 
flexibility  available  to  the  algorithm  than  in  the  one-way 
design  procedure,  and  it  generally  performs  more  poorly 
than  the  one-way  algorithm,  except  in  a  limited  range  of 
goal  psf  widths.  The  goal  psf  width  for  the  two-way  simu¬ 
lation  presented  in  Fig.  4  falls  within  this  range. 

Fig.  5  depicts  the  results  obtained  when  the  MSSE  al¬ 
gorithm  was  implemented  in  one-way  broadband  simula¬ 
tions.  It  can  be  seen  from  Figs.  5(a)  and  5(b)  that  the  goal 
and  system  psfs  are  qualitatively  quite  similar.  The  rela¬ 
tive  rms  error  was  0.003,  and  the  relative  peak  error  was 
0.05. 

Results  from  the  two-way  broadband  simulation  are 
shown  in  Fig.  6.  The  relative  rms  error  was  0.001,  and 
the  relative  peak  error  was  0.022. 

The  results  obtained  in  the  one-  and  two-way  broad¬ 
band  design  simulations  (Figs.  5  and  6)  must  be  inter¬ 
preted  with  caution  because  the  goal  psfs  used  are  difficult 
to  realize  using  spherical  waves,  as  can  be  seen  in  Figs.  5(a) 
and  6(a).  The  goal  psfs  had  flat  wavefronts  and  would  have 
been  easy  to  generate  using  plane  waves,  but  plane  waves 
are  an  unrealistic  model  of  the  ultrasonic  field  emitted  by 
transducer  elements.  Medical  ultrasonic  imaging  is  typi¬ 
cally  performed  in  the  near  field  of  the  transducer.  Because 
of  this,  both  transmitted  and  received  wavefronts  should 
properly  be  considered  as  spherical  wavefronts.  In  spite  of 
the  challenging  goal  psfe  used  here,  there  is  a  very  good 
qualitative  agreement  between  the  goal  and  system  psfe. 
Unlike  in  the  CW  simulations,  the  relative  errors  in  the 
two-way  broadband  simulation  are  smaller  than  that  ob¬ 
tained  in  the  one-way  simulation.  We  believe  that  this  is 
due  to  the  more  realistic  broadband  simulations  performed 
using  Field  II.  Field  II  uses  elements  of  finite  spatial  extent, 
as  opposed  to  the  ideal  point  sources  used  in  the  CW  sim¬ 
ulations.  The  use  of  such  elements  simulates  the  effects  of 
nonuniform  element  angular  response,  which  we  hypothe¬ 
size  to  be  partly  responsible  for  the  lower  two-way  relative 
errors.  The  size  of  the  window  of  analysis  and  the  choice 
of  spatial  and  temporal  sampling  rates  also  may  play  a 
significant  role.  However,  in  both  one-way  and  two-way 
simulations,  wc  were  able  to  approximate  a  very  challeng¬ 
ing  goal  psf. 

The  ability  of  the  MSSE  technique  to  improve  the  DOF 
is  illustrated  in  Fig.  7.  The  DOF  was  defined  in  terms  of 
correlation  coefficient  as  the  axial  region  over  which  the  co¬ 
efficient  remained  above  0.99.  In  the  CW  case,  the  DOF  in 
control  simulations  was  7.8  mm,  and  it  increased  by  249% 
to  27.2  mm  when  the  MSSE  technique  was  applied.  In  the 
broadband  simulations,  the  DOF  increased  from  4.1  mm 
for  control  to  17.4  mm  upon  the  application  of  the  MSSE 
technique,  an  increase  of  325%.  The  improvement  in  the 
DOF  in  CW  simulations  can  be  clearly  seen  in  Fig.  7(a), 
and  in  broadband  simulations  in  Fig.  7(c). 

The  DOF  also  was  defined  in  terms  of  the  FWHM  of 
the  mainlobe.  Here  we  considered  the  DOF  to  be  the  re¬ 
gion  within  which  the  FWHM  stayed  within  ±25%  of  its 
value  at  the  focus.  The  DOF  calculated  using  the  FWHM 
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criterion  increased  from  11  mm  in  control  simulations  to 
almost  the  entire  range  of  interrogation,  i.e.,  50  mm.  This 
represents  a  355%  increase  in  the  DOF  when  the  MSSE 
technique  was  applied.  The  CW  control  and  MSSE  tech¬ 
nique  FWHM  results  can  be  seen  in  Fig.  7(b).  Fig.  7(d) 
displays  the  FWHM  results  for  broadband  simulations.  In 
the  broadband  case,  the  DOF  evaluated  using  the  FWHM 
increased  from  12  mm  to  about  25.4  mm  upon  applica¬ 
tion  of  the  MSSE  technique,  an  increase  of  about  112%.  It 
can  be  seen  clearly  in  Figs.  7(b)  and  7(d)  that  the  main- 
lobe  FWHM  at  the  focus  was  identical  for  the  control  and 
MSSE  technique  cases,  but  it  varied  more  slowly  away  from 
the  focus  when  the  MSSE  method  was  applied.  Therefore, 
a  significant  improvement  in  the  DOF  was  obtained  in 
both  CW  and  broadband  simulations. 

A  more  qualitative  assessment  of  the  efficacy  of  the 
MSSE  technique  can  be  made  using  Figs.  7(c)  and  7(f). 
Fig.  7(e)  shows  the  CW  psfs  obtained  in  control  simula¬ 
tions  and  simulations  that  implemented  the  MSSE  tech¬ 
nique.  The  psfs  are  displayed  as  a  function  of  range  and 
lateral  position.  It  can  be  seen  that  there  is  significant 
broadening  of  the  psf  mainlobe  with  range  in  the  control 
simulation.  Fig.  7(f)  clearly  demonstrates  the  dramatic  im¬ 
provement  in  the  DOF  obtained  using  the  MSSE  beam¬ 
forming  technique.  This  improvement,  though,  comes  at 
the  cost  of  slightly  increased  sidelobe  levels.  This  is  due 
to  the  fact  that  the  MSSE  algorithm  minimizes  the  en¬ 
ergy  in  the  difference  between  the  goal  and  system  psfs. 
Because  most  of  the  energy  is  contained  in  the  mainlobe, 
the  algorithm  preferentially  optimizes  mainlobe  width  at 
the  cost  of  higher  sidelobe  levels.  However,  as  described  in 
[1],  we  can  apply  a  weighting  function  to  selectively  em¬ 
phasize  sidelobes  during  the  MSSE  design  process.  The 
procedure  involves  the  selection  of  an  appropriate  func¬ 
tion  that  weights  the  sidelobes  of  the  psf  more  than  the 
mainlobe  during  the  SSE  minimization  operation.  Imple¬ 
menting  the  MSSE  technique  then  would  yield  aperture 
weightings  that  preferentially  maintain  low  sidelobe  lev¬ 
els,  at  the  expense  of  mainlobe  width.  This  would  enable 
a  compromise  between  the  overall  performance  of  the  al¬ 
gorithm  and  sidelobe  levels.  However,  as  can  be  seen  from 
Fig.  7,  the  MSSE  technique  performs  much  better  overall 
than  the  conventional  techniques  of  dynamic  apodization 
and  dynamic  receive  focusing,  and  greatly  enhances  the 
depth  of  field. 

Fig.  8  illustrates  the  use  of  the  MSSE  design  method 
in  the  TAA.  Fig.  8(a)  shows  the  correlation  coefficients 
obtained  by  correlating  the  nontranslated  geometry  psf 
with  the  translated  geometry  psf  in  CW  simulations.  As 
previously  described,  the  correlation  DOF  was  defined  as 
the  axial  region  within  which  the  correlation  coefficient 
remained  over  0.99.  The  correlation  DOF  increased  from 
5.2  mm  in  control  simulations  to  46.9  mm  upon  the  ap¬ 
plication  of  the  MSSE  technique.  This  represents  an  802% 
increase  in  the  DOF.  Fig.  8(b)  displays  the  correlation  co¬ 
efficients  obtained  in  broadband  simulations.  The  DOF  in¬ 
creased  from  1.0  mm  to  8.8  mm  in  broadband  simulations, 
an  increase  of  780%. 


The  goal,  control,  and  MSSE  technique  psfs  at  ax¬ 
ial  ranges  of  4  mm,  6.5  mm,  and  9  mm  are  depicted  in 
Figs.  8(c),  8(d),  and  8(e),  respectively.  It  can  be  seen 
clearly  that  the  MSSE  technique  preferentially  optimizes 
the  mainlobe  of  the  obtained  psf,  at  the  cost  of  higher  side- 
lobes.  If  sidelobe  levels  are  too  high,  the  weighted  MSSE 
technique  described  in  [1]  can  be  used  to  achieve  a  com¬ 
promise  between  optimizing  the  mainlobe  and  controlling 
sidelobe  levels. 

The  MSSE  design  method  worked  exceedingly  well  over 
the  range  of  conditions  considered  in  this  paper,  but  we 
must  exercise  some  caution  in  interpreting  these  results. 
The  results  shown  in  Fig.  2  clearly  demonstrate  that  use 
of  an  appropriately  large  window  of  analysis  is  critical.  An¬ 
other  concern  is  the  effect  of  assuming  an  incorrect  prop¬ 
agation  model  in  the  derivation  of  the  optimum  aperture 
weights.  Errors  such  as  a  mismatch  in  the  assumed  and 
actual  wave  propagation  speeds  have  been  shown  to  have 
an  adverse  effect  on  the  design  method,  although  observed 
errors  weTe  small.  Phase  aberration  will  also  adversely  im¬ 
pact  the  performance  of  the  MSSE  algorithm,  because  it  is 
unaccounted  for  in  the  propagation  models  used.  Blocked 
or  dead  elements  will  also  have  an  effect,  because  any  as¬ 
sumed  propagation  functions  for  these  elements  will  cause 
an  undesired  contribution  to  the  ultrasound  field.  These 
effects  remain  to  be  investigated  in  detail,  but  the  initial 
simulation  results  that  have  been  presented  suggest  that 
the  MSSE  technique  is  a  robust  beamformer  design  tool. 

The  calculated  weights  do  not  lend  themselves  to  intu¬ 
ition,  but  the  MSSE  algorithm  may  be  considered  to  be 
analogous  to  the  discrete  Fourier  transform  [9]  or  wavelet 
transform  ([9,  pp.  500-502]).  The  discrete  Fourier  trans¬ 
form  expresses  a  signal  in  terms  of  multiple  narrowband 
signals  with  unique  frequencies.  We  can  process  (weight) 
these  narrowband  signals,  then  use  the  inverse  discrete 
Fourier  transform  to  sum  them  and  construct  a  new  signal 
that  we  desire.  Similarly,  the  wavelet  transform  decom¬ 
poses  a  signal  into  a  set  of  weighted  wavelets  that  can  be 
summed  to  construct  a  new  desired  signal.  The  MSSE  al¬ 
gorithm  is  analogous  because  it  decomposes  the  system 
response  into  the  contributions  of  the  individual  elements. 
These  individual  element  responses  then  are  weighted  and 
summed  to  construct  the  desired  system  response.  The  in¬ 
dividual  element  responses  are  not  necessarily  orthogonal 
such  as  the  kernel  functions  used  in  Fourier  and  wavelet 
decompositions,  but  the  MSSE  algorithm  operates  under 
the  same  principle. 

Simulations  have  shown  the  MSSE  beamformer  design 
technique  to  be  able  to  design  apertures  for  applications 
that  require  arbitrary  system  responses.  One  such  applica¬ 
tion  is  multidimensional  blood  velocity  estimation  as  de¬ 
scribed  in  [14]  and  [15],  The  specialized  psfs  required  for 
the  methods  in  [14]  and  [15]  can  be  optimally  generated 
using  the  MSSE  technique.  Because  the  MSSE  technique 
derives  aperture  weights  that  minimize  the  SSE  between 
the  goal  and  system  psfs,  the  calculated  weights  generate 
a  system  response  that  is  optimally  similar  to  the  goal 
response.  This  eliminates  the  usual  need  for  iteration  to 
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obtain  an  adequate  response.  The  technique  has  also  been 
shown  to  be  successful  in  developing  apertures  for  com¬ 
mon  beamformer  design  problems  such  as  limited  depth 
of  field.  Overall,  the  MSSE  technique  has  the  potential  to 
improve  beamforming  in  general,  with  much  better  control 
of  beam  parameters  than  is  possible  with  current  beam¬ 
forming  techniques.  The  direct  solutions  provided  by  this 
approach  also  have  the  potential  to  save  a  great  deal  of 
time  by  obviating  iterative  design. 

IV.  Conclusions 

The  MSSE  technique  has  been  shown  to  be  effective  in 
designing  ultrasound  systems  that  generate  arbitrary  de¬ 
sired  system  responses.  Simulation  results  in  one-way  and 
two-way  CW  and  broadband  systems  demonstrate  that 
it  is  straightforward  to  implement  and  can  be  applied  to 
a  wide  range  of  potential  applications.  Simulation  results 
obtained  by  implementing  the  beamforming  method  in  ex¬ 
amples  of  application  demonstrate  the  success  of  the  tech¬ 
nique  in  solving  common  problems  that  are  encountered 
in  ultrasound  imaging,  such  as  a  restricted  depth  of  field. 

Therefore,  the  MSSE  technique  has  been  shown  to  have 
significant  potential  to  improve  ultrasound  beamforming 
and  can  be  applied  in  any  ultrasound  application  in  which 
better  control  of  beam  parameters  is  desired.  Specifically, 
applications  that  require  specialized  psfs  are  well  suited 
for  the  technique.  There  is  no  iteration  involved;  there¬ 
fore,  design  time  is  considerably  reduced.  Further  inves¬ 
tigation  is  required  to  examine  the  effects  of  phase  aber¬ 
ration,  blocked  elements,  and  imposing  constraints  on  the 
calculated  aperture  weights.  But  our  simulations  indicate 
that  the  MSSE  technique  consistently  outperforms  current 
beamforming  techniques. 
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Abstract  -  Traditional  ultrasound  systems  measure 
backscatter  in  B-mode,  capturing  only  the  acoustic 
energy  that  is  reflected  directly  from  the  target 
region  to  the  transducer  face.  These  systems  fail  to 
utilize  the  information  in  the  echo  field  that  is 
scattered  in  other  directions  and  therefore  cannot 
characterize  the  angular  scattering  behavior  of  the 
targets  being  observed.  Since  target-specific 
angular  scattering  has  great  potential  as  a  source  of 
increased  contrast  in  biological  tissues,  it  is 
desirable  to  modify  the  method  of  acquisition  in 
order  to  obtain  reliable  information  about  this 
behavior.  However,  prior  systems  used  to 
investigate  this  information  have  been  clinically 
unwieldy  and  statistically  inaccurate  over  small 
regions.  We  have  implemented  a  method  of 
acquisition-  that  utilizes  the  translating  apertures 
algorithm  (TAA)  to  reliably  separate  target-specific 
angular  scatter  information  from  the  effects  of 
changing  acquisition  geometry.  This  acquisition 
method  has  been  implemented  in  real-time  on  a 
clinical  linear  array  system.  Seven  interrogation 
angles  are  acquired  for  each  imaging  line,  and  the 
TAA  is  implemented  repeatedly  across  the  array  to 
yield  per-pixel  maps  of  angular  scatter  behavior. 
We  present  comparisons  of  per-pixel  angular 
scatter  behavior  for  a  variety  of  target  types, 
including  correlation  analysis  of  a  Rayleigh-regime 
wire  target  phantom  and  comparative  image 
analysis  with  phantoms  containing  targets  of 
varying  compressibility  and  density.  It  is  shown 
that  the  TAA  maintains  a  high  per-pixel  correlation 
level  over  a  broad  range  of  interrogation  angles. 
Comparative  angular  scatter  imaging  is  shown  to 
yield  relative  contrast  improvements  on  the  order  of 
10-1 5dB  in  some  targets. 

I.  INTRODUCTION 

Traditional  clinical  ultrasound  imaging  methods 
yield  information  about  the  echogenicity  of  targets 


within  the  imaging  field  by  processing  only  those 
echoes  which  are  returned  directly  to  the  point  Of 
transmission  (backscatter).  Such  systems  do  not 
consider  portions  of  the  returning  echo  field  which 
are  scattered  in  other  directions,  and  thus  cannot 
provide  information  about  the  angular  scatter 
profile  of  insonified  targets.  The  character  of  the 
angular  scatter  profile  has  been  shown  to  contain 
significant  information  about  the  type  of  target 
being  observed  [1],  so  it  is  desirable  to  develop  an 
imaging  system  which  can  process  this  type  of 
information  in  a  clinically  useful  manner. 

In  general,  the  amplitude  of  the  echo  field  emitted 
by  a  Rayleigh  scattering  target  will  exhibit  an 
angular  dependence  that  is  proportional  to  its 
background-relative  compressibility  and  density: 


3  rc  L  k  2p,-p 


(1) 


where  a  is  the  scatterer  radius,  r  is  the  target  range, 
Kt  and  k  are  the  target  and  background 
compressibilities,  respectively,  pt  and  p  are  the 
target  and  background  densities,  respectively,  and 
Gs  is  the  scattering  angle  relative  to  backscatter 
(180°).  Generated  echoes  consist  of  the  summed 
contributions  of  an  omni-directional  wave  caused 
by  local  compressibility  variations  and  an  angle- 
dependent  (dipolar)  wave  caused  by  local  variations 
in  density,  but  this  information  cannot  be  separated 
with  a  conventional  imaging  geometry. 

We  have  developed  an  imaging  system  using  the 
translating  apertures  algorithm  (TAA)  [2]  which 
allows  for  the  reliable  evaluation  of  the  echo  field 
at  multiple  interrogation  angles.  This  allows  for  the 
omni-directional  and  angle-dependent  components 
of  the  echo  field  to  be  evaluated  separately,  which 
introduces  compressibility  and  density  variations  as 
a  new  source  of  potential  image  contrast.  This 
imaging  method  is  implemented  in  real-time  on  a 
GE  Logiq  700MR  clinical  imaging  system  using  a 
7.5MHz  linear  array  probe.  We  present  the  initial 
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evaluation  of  this  system’s  ability  to  isolate  angular 
scatter  information  through  analysis  of  a  variety  of 
target  types. 

II.  METHODS 

The  advantage  of  the  translating  apertures 
algorithm  over  traditional  angular  scatter 
acquisition  geometries  [3]  is  the  reliable  isolation  of 
target-specific  echo  information  from  other  system 
effects.  Moving  the  transmit  and  receive  apertures 
in  equal  and  opposite  directions  to  increase 
interrogation  angle  provides  a  stable  system  point 
spread  function  through  a  broad  angular  range. 

In  order  to  implement  the  TAA  on  a  linear  array  in 
a  useful  manner  that  can  generate  images,  it  is 
necessary'  to  modify  system  behavior  during 
transmission  and  reception  such  that  a  variety  of 
angles  can  be  interrogated  for  every  point  in  lateral 
space  along  the  array.  Implementation  of  the  TAA 
on  the  Logiq.  700  system  involves  extensive 
modification  to  allow  for  precise  control  of  the 
apertures  being  utilized  for  transmission  and 
reception  for  every  set  of  pulses  that  the  system 
fires.  Since  linear  array  systems  are  designed  to 
focus  at  multiple  depths  and  (typically)  fire  one  set 
of  focused  pulses  per  lateral  image  line  per  focus, 
modifications  must  be  made  to  interrogate  a  single 
region  in  space  at  multiple  angles. 

Firstly,  identical  transmit/receive  beamforming  is 
implemented  across  multiple  system  focal  zones, 
such  that  the  system  interrogates  the  same  spatial 
region  many  times  (according  to  the  number  of 
interrogation  angles  desired).  It  should  be 
emphasized  that  though  the  system  interprets  each 
of  these  redundant  firings  as  a  different  spatial 
focal  zone,  all  the  calculated  time  delays  are 
identical  for  each  zone.  Additionally,  dynamic 
receive  focusing  and  transmit/receive  apodization 
(which  are  depth  dependent)  are  disabled.  Aperture 
translation  is  achieved  through  pre-calculated 
electronic  channel  maps  which  can  independently 
turn  any  array  element  on  or  off  during  transmit  and 
receive  operations.  By  using  the  same 
beamforming  calculations  for  every,  “zone,”  small 
angles  can  be  interrogated  by  turning  on 
transmit/receive  elements  near  the  center  of  the 
active  aperture,  and  large  angles  can  be  interrogated 
by  enabling  elements  near  the  edge.  All  imaging 
vectors  are  placed  in  line  with  the  physical  imaging 


elements  to  assure  absolute  symmetry  as  TAA  is 
applied  to  separate  the  active  transmit  and  receive 
apertures  in  equal  and  opposite  directions.  The 
firing  order  of  the  system  for  a  four-angle 
acquisition  thus  looks  like  this: 
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Figure  1 :  Multi-angle  acquisition  on  a  linear  array 

This  example  demonstrates  acquisition  of  a 
backscatter  angle  (0,)  and  three  other  interrogation 
angles  for  two  different  imaging  lines  (A  and  B). 
The  firing  order  for  the  system  would  be  A0i,  A02, 
A03,  A04,  B0j,  B02,  B03,  B04.  There  are 
approximately  200  lateral  spatial  vectors  per  image 
frame,  and  up  to  seven  interrogation  angles  can  be 
interrogated  per  spatial  imaging  line. . 

Complex  echo  information  is  obtained  in  the  form 
of  summed  IQ  data  that  is  offloaded  from  the  Logiq 
700  to  a  pc-based  storage  system,  where  it  is  then 
unpacked/de-interlaced  and  processed  for  the 
purposes  of  comparative  evaluation. 

III.  RESULTS 

Rayleigh-regime  wire  target 
As  it  is  essential  that  the  TAA  acquisition 
maintain  a  highly  correlated  system  point  spread 
function  across  all  angles  of  interrogation,  an  initial 
experiment  was  performed  using  a  20pm  diameter 
stainless  steel  wire(<  X/IO,  approximating  a 
Rayleigh  scattering  target)  in  a  degassed  water 
bath.  Seven  acquisition  angles  were  acquired,  from 
backscatter  (180°)  to  162°  in  3-degree  increments. 
Since  transmit  and  receive  apodization  could  not  be 
controlled  independently,  no  apodization  was  used 
for  any  of  these  experiments.  Several  parameters 
were  varied  to  test  their  effect  on  psf  correlation 
levels,  including  the  size  of  the  transmit/receive 
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apertures,  the  speed  of  sound  used  by  the  system  to 
calculate  time  delays,  and  the  distance  the  wire 
target  was  placed  away  from  the  transducer  face. 
In  general,  varying  the  calculated  speed  of  sound 
produced  the  same  results  as  moving  the  transducer 
axially  toward/away  from  the  wire  target,  and 
larger  transmit/receive  apertures  provided  better 
correlation  values  (largely  due  to  increased  SNR). 
Using  an  8  element  (1.6mm)  aperture  with  the 
system  focused  at  12  mm,  the  following  plot  was 
generated  to  demonstrate  correlation  depth  of  field: 


Figure  2.  Correlation  coefficient  vs.  distance  from 
transducer  face 

The  complex  correlation  coefficient  p  was 
calculated  relative  to  the  acquisition  at  180°  (note 
that  the  angle  0  on  the  legend  corresponds  to  1 80- 
Os,  it  is  more  intuitive  to  refer  to  this  angle  since 
increasing  0  indicates  an  increased  separation 
between  transmit  and  receive  apertures)  .  Looking 
at  the  absolute  value  of  the  complex  correlation 
coefficient  it  is  clear  that  correlation  drops  when 
interrogation  angle  is  increased,  as  well  as  when  the 
transducer  is  moved  farther  from  the  focus  at 
12mm.  All  angles  exhibit  high  correlation  levels 
at/near  the  focus,  and  correlation  remains  above 
0.95  for  approximately  5  mm  for  all  angles  of 
interrogation. 

The  simplest  method  of  comparative  imaging  that 
can  be  employed  to  quantify  differences  in  angular 
scatter  is  common-  and  difference-  (c-  and  d-) 
weighted  subtractive  imaging  [2].  A  difference- 
weighted  image  is  formed  via  the  subtraction  of  the 
complex  echoes  acquired  at  an  angle  of  interest 
from  those  acquired  at  the  backscatter  angle.  This 
eliminates  the  portion  of  the  echo  field  common  to 
every  acquisition  angle  and  emphasizes  those 
portions  of  the  echo  field  which  change 


significantly  as  the  interrogation  angle  is  increased. 
Subtracting  the  difference  echoes  from  the  echoes 
acquired  at  the  separation  angle  of  interest  yields 
the  common-weighted  echo  set.  Common- 
weighted  images  emphasize  targets  dominated  by 
local  variations  in  compressibility,  while 
difference-weighted  images  emphasize  local 
variations  in  density. 

To  evaluate  the  potential  efficacy  of  a 
straightforward  subtraction  image  on  the  system, 
the  normalized  difference  energies  (DE)  of  all  data 
sets  were  calculated  as  follows: 


DE  —  A0  /  Aq  A 0  —  Aq  /  ^ Aq  Aq 


(2) 


Where  Ao  is  the  complex  backscatter  echo  data  and 
Ae  is  the  complex  echo  data  a  higher  angle  of 
interrogation.  Initial  evaluation  of  this  metric  is 
shown  in  Fig.  3: 


Figure  3.  Normalized  difference  energy  vs. 
distance  from  focus 


Results  showed  a  higher  than  expected  rise  in 
difference  energy  at  high  interrogation  angles,  so  in 
order  to  improve  the  depth  of  field  over  which  C- 
and  D-  imaging  could  be  successfully  employed, 
receive  time  warping  (RTW)  was  applied  to  align 
the  data  acquired  at  backscatter  to  the  echoes 
acquired  at  other  angles.  RTW  compensates  for  the 
relative  difference  in  pathlength  between  the  focus 
and  other  points  in  the  imaging  field  when 
considering  angular  acquisition  vs.  backscatter.  A 
complex  phase  rotation  is  applied  to  the  angular 
acquisition  data  sets  to  make  up  for  this  difference 
(Fig.  4).  RTW  greatly  enhanced  the  effective  depth 
of  field  for  C-  and  D-  type  image  comparisons, 
offering  an  order  of.  magnitude  improvement  at 
larger  interrogation  angles. 
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Figure  4.  Normalized  difference  energy  vs.  distance 
from  transducer  face,  receive  time  warping  active 


To  quantify  the  effectiveness  of  C-  and  D- 
imaging  on  different  target  types,  a  phantom  was 
constructed  featuring  three  different  types  of 
200pm  diameter  wire.  The  base  of  the  phantom 
was  high-concentration  gelatin  with  50pm  diameter 
Sephadex  suspended  within  to  provide  background 
scattering.  This  gelatin  was  poured  upon  a  steel 
wire,  nylon  monofilament,  and  cotton/polyester 
thread  placed  approximately  1  cm  apart,  and  were 
set  to  bisect  the  imaging  plane  of  the  transducer  at 
the  focal  depth  of  12  mm.  B-Mode  (at  0  =  0°)  C-, 
and  D-weighted  images  (from  0  =  9.679  degrees) 
are  presented  in  Fig.  5. 


B-Mode  Image 


Figure  5.  B-Mode,  C-weighted,  and  D-weighted 
image  of  steel  (left),  nylon  (center),  and  cotton 
(right)  wire  targets 

Average  pixel  brightness  in  the  three  wire  regions 
was  compared  to  that  of  a  neutral  background 
region  also  at  the  focus.  In  the  B-Mode  images  the 


contrast  between  the  wires  and  the  background  is 
nearly  the  same:  21.3dB,  21.9dB,  and  23.8dB 
respectively  for  the  steel,  nylon,  and  cotton  wires. 
Note  that  they  are  difficult  to  differentiate  from  the 
B-Mode  image  alone.  In  the  C-weighted  image  we 
see  contrast  levels  of  19.4dB,  20.5dB,and  25.3dB 
(note  the  enhanced  contrast  for  the  highly 
compressible  cotton  thread).  In  the  D-weighted 
image  we  see  contrast  levels  of  27.90dB,  22.50dB, 
and  17.07dB  for  the  steel,  nylon,  and  cotton 
threads,  such  that  all  three  wire  types  are  clearly 
differentiable  in  the  image. 


IV.  CONCLUSION 

Although  there  are  a  few  considerable  limitations 
currently  involved  with  implementing  TAA  on  a 
linear  array  system  (most  notably,  the  inability  to 
independently  apodize  the  transmit  and  receive 
apertures),  it  has  been  shown  that  there  is  great 
potential  in  this  technique  for  providing  new 
sources  of  image  contrast  from  previously 
unobservable  target  characteristics.  Further 
refinement  of  the  technique  should  allow  for  the 
generation  of  higher  quality  overall  image  data,  as 
well  as  the  potential  real-time  displa>  of 
comparative  (C/D  -type)  images. 
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Abstract  -  Current  aperture  design  techniques  do  not 
allow  for  the  design  of  apertures  that  produce  a  beam 
pattern  optimally  similar  to  the  desired  pattern.  A 
flexible  beamforming  technique  that  enables  the 
optimal  design  of  apertures  for  a  desired  system 
response  is  presented.  This  technique  involves  a 
linear  algebra  formulation  of  the  Sum  Squared  Error 
(SSE)  between  the  point  spread  function  (psf)  of  the 
system,  and  an  ideal  or  desired  psf.  Minimization  of 
this  SSE  yields  the  optimum  aperture  weightings.  A 
brief  overview  of  the  application  of  the  technique  for 
some  common  design  objectives,  along  with 
simulation  results  is  also  presented. 

I.  INTRODUCTION 

A  common  task  in  ultrasound  imaging  is  the 
design  of  apertures  for  either  a  specific  imaging 
application,  or  to  improve  performance  in  an  existing 
application.  This  is  usually  performed  in  an  iterative 
and  ad-hoc  manner.  These  techniques  may  result  in  a 
system  response  close  to  the  desired  response,  but 
they  do  mot  guarantee  optimization  of  ultrasound 
beam  parameters  such  as  main-lobe  width  and  side- 
lobe  levels.  Also,  given  a  desired  beam  pattern,  it  is 
not  possible  to  design  apertures  that  produce  a  beam 
pattern  that  optimally  resembles  the  desired  beam. 

We  propose  a  general  aperture  design  method, 
with  rigorous  theory,  that  can  be  applied  in  arbitrary 
system  geometries  to  design  apertures  optimizing 
beam  parameters.  Our  technique  involves  a  linear 
algebra  formulation  of  the  Sum  Squared  Error  (SSE) 
between  the  system  point  spread  function  (psf)  and 
the  desired  or  ideal  psf.  Minimization  of  this  error 
yields  .unique  aperture  weightings  that  force  the 
system  psf  to  resemble  the  desired  psf.  It  is  similar  to 
the  technique  used  by  Ebbini  [1],  and  by  Ebbini  and 
O’Donnell  [2].  Our  Minimized  Sum  Squared  Error 
(MSSE)  technique  is  more  general,  however,  because 
it  enables  the  use  of  arbitrary  propagation  functions. 


Another  distinction  is  that  in  [1]  and  [2],  a  few 
control  points  were  used  in  order  to  ensure  an 
underdetermined  system  of  equations  and  obtain  an 
exact  beam  pattern  at  those  few  points,  while  we  use 
the  entire  ideal  psf  to  obtain  the  least  squares  solution 
of  an  overdetermined  system  of  equations.  This 
method  enables  excellent  control  of  the  system  psf, 
and  has  a  significant  impact  on  aperture  design  for 
several  applications  such  as  improved  depth  of  field. 

II.  THEORY 

The  phase  and  amplitude  of  the  ultrasonic  field  at 
a  point  in  space  due  to  a  transducer  element  depends 
on  several  factors.  These  include  the  Euclidean 
distance  between  the  point  and  the  element,  the 
orientation  of  the  element  relative  to  the  point,  the 
frequency  of  the  emitted  wave,  and  frequency 
dependent  attenuation  of  the  medium,  assuming  linear 
propagation.  The  complex  field  at  the  point  under 
consideration  can  be  written  as  the  product  of  a 
propagation  function,  s ,  which  incorporates  any  or 
all  of  the  above  mentioned  factors,  and  the  weighting 
(possibly  complex)  applied  to  the  element,  w.  i.e. 
sw .  The  one-way  M-point  lateral  point  spread 
function  (psf)  at  the  range  z  can  be  represented  as, 


where  Sz  is  an  M  x  N  matrix  of  propagation 
functions  in  which  each  element  su  is  the 
propagation  function  that  determines  the  field  at  a 
point  i  due  to  element  j ,  W  is  an  N  x  1  vector  of 

aperture  weightings  in  which  each  element  Wj  is  the 
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weighting  applied  to  the  j th  element,  and  Pz  is  the 
resulting  psf,  which  is  an  M  x  1  vector.  This 
formulation  permits  analysis  with  complicated, 
propagation  functions  that  may  include  limited 
element  angular  response,  frequency  dependent 
attenuation,  and  other  difficult  to  model  factors. 

Using  equation  1  and  by  applying  the  well-known 
Radar  Equation,  the  two-way  psf  can  be  written  as, 

PTK=PnP*,=(S!T)(SlR)  (2) 

where  Pyz  and  P^  are  the  one-way  transmit  and 

receive  psfs  respectively,  T  and  R  are  the  transmit 
and  receive  aperture  weightings,  and  V  indicates 
point  by  point  multiplication.  The  propagation 
function  is  the  same  on  transmit  and  receive  due  to 
acoustic  reciprocity.  Equation  2  can  be  rewritten  as, 

Pm  =  Pfzd  Prz  =  Pnd  R  =  PTldS  R ,  (3) 


receive  aperture  to  obtain  a  system  psf,  PTRi ,  at  the 
range  z ,  that  optimally  resembles  the  desired  or  ideal 
psf,  PTRji ..This  aperture  design  method  guarantees 
optimal  beam  patterns.  We  describe  some  common 
design  objectives,  and  the  application  of  our  method 
to  design  apertures  that  achieve  these  objectives,  to 
demonstrate  the  effectiveness  of  the  MSSE  technique. 

Objective:  Enhanced  Depth  of  field  -  The  depth  of 
field  (DOF)  of  an  ultrasound  imaging  system  is 
generally  defined  as  the  axial  region  over  which  the 
system  is  in  focus,  or  more  rigorously,  the  axial 
region  over  which  the  system  response  is  uniform 
within  some  predetermined  limit.  It  is  generally 
desired  that  the  system  psf  remains  similar  to  the  psf 
at  the  focus  for  as  large  an  axial  range  as  possible. 
Current  state  of  the  art  techniques  to  improve  depth  of 
field  include  transmit  apodization,  dynamic  receive 
apodization  and  dynamic  receive  focusing. 


where  PTzd  is  a  diagonal  M  x  M  matrix  with  the 
elements  of  PTl  along  its  0Ul  diagonal,  and 
Phds  =  Pfzd  •  This  changes  the  point  multiplication 
operation  to  a  regular  matrix  multiplication  operation. 
We  can  characterize  the  degree  of  similarity  between 
the  psf  at  some  range  z ,  and  some  ideal  psf  by  the 
SSE  between  them.  Using  equation  3  the  SSE  is, 

SSE  =  [ptRi  -  Pmi  J  -  Pjjy  )  ^ 

=  {P'rzds  P  ~  Pr/tzi )  (frzds  R  ~  Pmi ) 

where  Pj-Rzi  is  the  ideal  psf,  and  the  superscript  t 
denotes  a  complex  conjugate  operation. 


This  is  formulation  is  well-known  in  signal 
processing,  and  using  [3],  we  can  obtain  a  set  of 
receive  weightings  to  be  applied  so  that  the  SSE 
between  the  generated  and  ideal  psfs  is  minimized. 


P  '  P 

rTidS  rTkii 


Tzds  1  mi 


(5) 


where  PTxdS*  is  the  pseudoinverse  of  PTldS  . 


III.  DISCUSSION 


Equation  5  specifies  the  complex  weightings  to 
be  applied  to  the  transducer  elements  constituting  the 


However  effective  the  above  techniques  are,  they 
are  ad-hoc  and  lack  formal  theory  describing  the 
effectiveness  in  improving  depth  of  field.  Our 
objective  is  to  derive  receive  weightings  that  force  the 
psf  at  each  specific  range  of  interrogation  to  be 
optimally  similar  to  the  psf  at  the  fociis,  and  use  these 
weightings  to  implement  dynamic  weighting  to 
maximize  depth  of  field.  This  can  be  done  easily  by 
setting  the  ideal  psf,  PTKjj ,  in  equation  5  to  be  the  psf 
at  the  focus.  We  will  obtain  receive  weightings  for 
the  range  of  interest,  z ,  that  will  generate  a  psf  that  is 
optimally  similar  to  the  psf  at  the  focus. 

Objective:  Correlation  Depth  of  Field  in  Translated 
Aperture  Geometries  —  We  have  proposed  using  the 
Translating  Apertures  Algorithm  (TAA)  as  the 
foundation  of  angular  scatter  imaging  methods  [4]. 
The  TAA  results  in  a  considerably  reduced  depth  of 
field  as  the  transmit  and  receive  apertures  are 
translated.  Our  technique  of  dynamic  receive- aperture 
weightings  can  be  applied  to  improve  the  correlation 
between  the  backscatter  (non-translated)  and  angular 
scatter  (translated)  psfs  at  each  range  of  interest,  and 
thereby  increase  the  correlation  depth  of  field.  We 
can  derive  the  receive  weightings  to  be  applied  in  the 
translated  apertures  geometry  to  maximize  the 
correlation  by  minimizing  the  SSE  between  the  two 
psfs.  We  can  write  the  two-way  psf  for  the  translated 
geometry  at  range  z  as  follows. 
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(9) 


where  Pmd  is  a  diagonalized  M  x  M  matrix  with  the 
elements  of  PTll  along  the  0*  diagonal,  the  subscript 
“1”  denotes  the  translated  geometry  and 
Pnids  =  Pmd  Sgzi  ■  The  SSE  between  the  psfs  for  the 
backscatter  and  angular  scatter  geometries  is, 


SSE  —  (Pffoi  P'nho  )  (.Ptrii  ^V/&o) 

=  (pTzldS  Pi  ~  PtKzoJ  (PnidS  P\  -Ptxzo) 


(10) 


corresponds  to  the  lateral  psf  at  a  single  axial  range. 
This  range  was  varied  from  0.31  cm  to  5  cm  and  was 
sampled  every  100 //m.  Figure  la  consists  of  the 
lateral  psfs  corresponding  to  the  control  case  with  no 
apodization.  Figure  lb  is  made  up  of  the  psfs 
obtained  when  dynamic  apodization  and  dynamic 
receive  focusing  were  applied,  along  with  a  range 
dependent  gain  function.  Figure  1c  consists  of  the 
psfs  obtained  when  our  MSSE  technique  was  applied. 
It  can  be  seen  that  beam  characteristics  were 
maintained  over  a  longer  range  for  our  technique, 
than  when  conventional  beamforming  techniques 
were  applied. 


The  receive  weightings  to  be  applied  Ult 
angular  scatter  geometry  that  yield  the  optimum 
correlation  depth  of  field  are  therefore  [3], 

^1  —  {^TzldS  P Tz\dS  )  ^TzUS  ^77tz0  =  ^Tz)dS 

where  PmdS*  is  the  pseudoinverse  of  PmdS 


The  above  design  objectives  illustrate  the 
flexibility  of  the  MSSE  technique.  The  method, 
however,  is  not  limited  to  these  examples  and  can  be 
used  to  design  apertures  to  obtain  any  arbitrary 
system  response. 

IV.  RESULTS  AND  DISCUSSION 


Simulations  of  the  design  method  were  performed 
by  implementing  code  in  Matlab.  We  used  a 
discretized  Rayleigh-Sommerfeld  formulation  to 
generate  the  propagation  matrices.  The  control 
parameters  are  described  in  table  1. 


Number  of  elements 

32 

Element  spacing 

200  microns 

Apodization 

Hann  window 

Focus 

1.2  cm 

Frequency  of  operation 

10  MHz 

Number  of  field  points 

351 

Field  point  spacing 

20  microns 

Table  1.  Control  parameters. 

Enhanced  Depth  of  Field  -  Figure  1  demonstrates  the 
effectiveness  of  our  technique  in  improving  the  depth 
of  field.  Every  column  in  each  of  the  three  images 


Correlation  Depth  of  Field  in  Translated  Aperture 
Geometries  -  Figure  2  illustrates  the  application  of 
the  MSSE  technique  in  translated  aperture 
geometries.  The  correlation  coefficients  were 
obtained  by  correlating  the  translated  aperture 
(shifted  by  10  elements)  psf  with  the  non-translated 
aperture  psf  at  the  same  range.  It  can  be  seen  that  the 
application  of  our  technique  results  in  a  significantly 
higher  correlation  than  the  control  case. 

Although  results  are  not  shown,  it  is  also  possible 
to  design  transmit  apertures  that  produce  limited 
diffraction  transmit  beams  and  maintain  transmit 
beam  characteristics  for  a  significantly  larger  range. 
This  can  be  done  by  tiling  the  one-way  psfs  at  the 
specific  axial  ranges  over  which  the  beam 
characteristics  are  to  be  maintained,  and  tiling  the 
ideal  psfs,  one  for  each  range  of  interest.  The  SSE 
between  the  tiled  ideal  and  actual  psfs  can  then  be 
formulated  using  equation  1.  Minimizing  this  SSE 
will  yield  the  transmit  weightings  that  produce  the 
optimized  limited  diffraction  transmit  beam. 

IV.  CONCLUSION 

The  Minimum  Sum  Squared  Error  (MSSE) 
technique  provides  a  general  method  for  the  design  of 
aperture  weightings  that  can  be  applied  in  arbitrary 
system  geometries  to  design  apertures  in  order  to 
obtain  a  beam  pattern  optimally  similar  to  the  desired 
pattern.  By  applying  the  technique  to  some  common 
design  challenges  in  medical  ultrasound,  it  has  been 
shown  to  have  the  potential  to  significantly  improve 
the  performance  of  imaging  systems. 
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Figure  1.  Images  obtained  from  lateral  psfs  at 
multiple  ranges.  Each  column  consists  of  the  psf  at  a 
single  range,  la  is  the  control  case,  lb  is  dynamic 
apodization  and  dynamic  receive  focusing  with  range 
dependent  gain,  and  lc  is  our  MSSE  technique. 


2  3 

Axial  Range  (cm) 


Figure  2.  Correlation  coefficients  obtained  by 
correlating  the  shifted  aperture  (10  element  shift)  psf 
with  the  non-shifted  aperture  psf  at  the  same  range. 

IV.  FUTURE  WORK 

The  technique  described  above  is  implemented 
using  a  continuous  wave  (CW)  formulation.  We  have 
also  adapted  the  technique  for  broadband  systems. 
We  are  currently  experimenting  with  broadband 
simulations,  and  limited  diffraction  transmit  beams. 
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ABSTRACT 

Conventional  techniques  used  to  design  transducer  apertures  for  medical  ultrasound  are  generally  iterative  and  ad-hoc. 
They  do  not  guarantee  optimization  of  parameters  such  as  mainlobe  width  and  sidelobe  levels.  We  propose  a  dynamic 
aperture  weighting  technique,  called  the  Minimum  Sum  Squared  Error  (MSSE)  technique,  that  can  be  applied  in 
arbitrary  system  geometries  to  design  apertures  optimizing  these  parameters.  The  MSSE  technique  utilizes  a  linear 
algebra  formulation  of  the  Sum  Squared  Error  (SSE)  between  the  point  spread  function  (psf)  of  the  system,  and  a  goal  or 
desired  psf.  We  have  developed  a  closed  form  expression  for  the  aperture  weightings  that  minimize  this  error  and 
optimize  the  psf  at  any  range.  We  present  analysis  for  Continuous  Wave  (CW)  and  broadband  systems,  and  present 
simulations  that  illustrate  the  flexibility  of  the  technique. 

Keywords:  ultrasound,  beamforming,  aperture  design,  minimum  sum  squared  error,  SSE,  MSSE 

1.  INTRODUCTION 

In  ultrasonic  imaging,  the  characteristics  of  the  ultrasound  beam  fundamentally  affect  the  quality  of  the  data  obtained, 
and  therefore  need  to  be  carefully  adjusted  to  obtain  the  desired  system  response.  Beam  parameters  such  as  the  mainlobe 
width  and  sidelobe  levels  can  be  adjusted  by  changing  the  amplitude  and  phase  (time  delay)  of  the  weightings  applied  to 
the  active  elements,  and  also  by  controlling  the  size  of  the  active  aperture  (the  number  of  active  elements)  and  the 
frequency  of  operation. 

These  beamforming  parameters  do  not  act  independently;  altering  one  changes  the  impact  of  each  of  the  others. 
Consequently,  beamformer  design  is  a  complicated  multiparameter  optimization  problem.  Because  of  this  complexity, 
beamformer  parameters  are  typically  determined  using  a  combination  of  ad  hoc  methods,  simplified  theory,  and  iterative 
simulation  and  experimentation.  While  these  methods  are  effective,  they  are  time  consuming  and  provide  no  guarantee 
that  an  optimized  solution  has  been  found. 

We  propose  a  general  aperture  design  method  that  can  be  applied  in  arbitrary  system  geometries  to  design  apertures  that 
optimize  beam  parameters.  Our  technique  utilizes  a  linear  algebra  formulation  of  the  Sum  Squared  Error  (SSE)  between 
the  system  point  spread  function  (psf)  and  the  desired  or  goal  psf.  Minimization  of  this  error  yields  unique  aperture 
weightings  that  force  the  system  psf  to  resemble  the  desired  psf.  It  is  similar  to  the  technique  used  by  Ebbini  et  al  [1]  to 
generate  specialized  beam  patterns  for  hyperthermia,  and  by  Li  et  al  [2]  for  the  compensation  of  blocked  elements.  There 
are  several  differences,  however,  between  these  methods  and  the  technique  we  describe.  The  analysis  in  [1]  and  [2]  uses 
only  a  few  control  points,  while  we  use  the  entire  system  psf  to  form  an  overdetermined  system  of  equations  that  we  then 
solve.  Another  important  distinction  is  that  unlike  [1]  and  [2],  we  present  analysis  for  both  CW  and  broadband  systems. 
This  paper  outlines  the  theoretical  description  of  the  MSSE  technique  for  narrowband  and  broadband  systems  and 
discusses  a  few  examples  of  application.  Simulation  results  for  these  examples  are  also  described. 

2.  THEORY 


2.1  One-way  Continuous  Wave  (CW)  formulation 

The  phase  and  amplitude  of  the  ultrasonic  field  at  a  point  in  space  due  to  an  ultrasound  transducer  element  depends  on 
several  factors.  These  include  the  Euclidean  distance  between  the  point  and  the  element,  the  orientation  of  the  element 
relative  to  the  point,  the  frequency  of  the  emitted  wave,  and  frequency  dependent  attenuation  of  the  medium.  Assuming 
linear  propagation,  the  one-way  M-point  lateral  point  spread  function  (psf)  at  the  range  z  can  be  represented  as  follows. 


where  S.  is  an  M  x  N  matrix  of  complex  propagation  functions  comprising  elements  of  the  form  Stj ,  which 
represents  the  propagation  function  that  determines  the  field  at  a  point  i  due  to  element  j .  W  is  an  N  x  1  vector  of 
aperture  weightings  in  which  each  element  Wj  is  the  weighting  applied  to  the  j th  element.  Pz  is  the  resulting  M  x  1 

psf.  This  formulation  permits  analysis  with  complicated  propagation  functions  that  may  include  limited  element  angular 
response  and  other  such  factors.  The  transmit  psf  at  the  range  z  can  therefore  be  expressed  as  follows, 

PTz  =S2T  (2) 

where  T  comprises  the  transmit  aperture  weightings.  Let  us  suppose  that  the  desired  one-way  system  transmit  psf  is  PTz 

for  the  application  of  interest.  We  can  then  characterize  the  degree  of  similarity  between  the  goal  psf  PTz  and  the  actual 
system  psf  PTz  by  the  Sum  Squared  Error  (SSE)  between  them. 

sse=(pT2-pTz)'  (pTz-PTz) 

=  (szt-pTz)'{szt-pTz)  (3) 

where  the  superscript  “  /  ”  denotes  a  conjugate  transpose  operation.  Minimizing  the  SSE  yields  the  transmit  weights  that 
produce  the  system  psf  that  is  optimally  similar  to  the  goal  psf.  The  formulation  in  equation  3  is  common  in  signal 
processing,  and  significant  literature  exists  on  the  solution  to  the  equation  with  the  minimum  SSE.  Using  [3]  the  transmit 
aperture  weightings  that  minimize  the  SSE  are  given  by, 


T  =  k'  szY'sz'pTz=sz»pTz 


where  S*  is  the  pseudoinverse  of  Sz .  Equation  4  describes  the  calculation  of  the  transmit  weightings  that  yield  the 
system  psf  at  the  range  z  that  is  optimally  similar  to  the  goal  psf. 

2.2  Two-way  Continuous  Wave  (CW)  formulation 

Using  the  analysis  in  the  previous  subsection  and  by  applying  the  well-known  RADAR  equation  [4],  the  two-way  psf  is, 

PTRz  =  PTz  '  PRz  ~  {$z  T)  (Sz  R),  (5) 

where  V  indicates  point  multiplication.  Equation  5  can  be  rewritten  as, 

PTRz  ~  PTzd  PRz  ~  PTzd  S.R  =  PTzds  R  ,  ((,) 


where  PTzJ  is  a  diagonal  M  x  M  matrix  with  the  elements  of  PTz  along  its  0,h  diagonal,  and  PTzds  =  PTzd  S2 .  If  PTR. 
is  the  goal  psf,  the  SSE  between  the  system  and  goal  psfs  can  be  expressed  in  a  similar  fashion  to  the  one-way 
formulation  in  equation  2.  The  receive  aperture  weights  that  minimize  the  SSE  can  then  be  determined  as  shown  below. 


Figure  1 .  Illustration  of  broadband  formulation.  The  pulse 
applied  to  each  element  is  convolved  with  a  set  of  weights, 
which  are  distinct  for  each  element.  This  operation  is 
analogous  to  implementing  an  FIR  filter  on  each 
element/channel. 


SSE  —  [PfRi  ~  ptrz  )  ^ 

[pTRz  -PtRz) 

(7) 

=  {PTzdS  R  ~  PTRz 

)  {PfzdS  P  ~  PtRz  ) 

R  “  (*7 'zdS  Pj zdS  )  Pitas  PtRz  ”  PlzdS  PtRz 

(8) 

Figure  2.  Illustration  of  the  exploitation  of  symmetry.  Since 
the  aperture  is  symmetric  about  the  center  elements,  pairs  of 
elements  that  would  have  the  same  weights  are  grouped 
together.  Also,  the  lateral  symmetry  of  the  psf  about  the 
center  axis  permits  analysis  with  just  one  half  of  the  psf. 


where  PTzdS*  is  the  pseudoinverse  of  PTzds .  Equation  8  specifies  the  complex  weightings  to  be  applied  to  the  transducer 
elements  constituting  the  receive  aperture  in  order  to  obtain  a  two-way  system  psf  PTRz  at  the  range  z  ,  that  optimally 
resembles  the  desired  or  goal  psf  ?TRz- 

2.3  One-way  broadband  formulation 

The  CW  formulation  described  above  will  have  limited  accuracy  in  the  analysis  of  broadband  systems.  For  this  case,  we 
have  developed  a  modified  formulation.  The  one-way  point  spread  function  (psf)  PTz  at  a  specific  range  z  is  a  function 
of  lateral  position  and  time,  and  can  be  represented  as  follows. 

PTz  =AzT  (9) 

where  Az  is  a  propagation  matrix  that  depends  on  the  excitation  pulse  and  the  impulse  responses  of  the  elements 
comprising  the  transmit  aperture.  It  is  a  function  of  time  and  the  spatial  positions  of  the  element  and  field  point  under 
consideration.  It  describes  the  contribution  of  each  element  at  each  field  point  as  a  function  of  time.  The  generation  of  a 
one-way  psf  is  shown  in  figure  1 .  The  SSE  between  the  system  psf  and  the  desired  psf  can  be  expressed  as, 


(10) 


sse  =  (pTz-pT:)t{pTz-pTz), 

where  the  superscript  “  T  ”  denotes  a  transpose  operation.  Using  equation  9,  we  can  solve  for  the  transmit  weightings 
that  minimize  the  SSE  in  equation  10. 

t=(aJ  Azy  AjpTi=Az»pTz  (11) 

2.4  Two-way  broadband  formulation 

Similar  to  the  one-way  psf,  the  two-way  pulse-echo  psf  can  also  be  expressed  in  a  linear  algebra  formulation. 

PTRz=AzzR> 


where  A s  is  a  function  of  the  transmit  aperture  weights,  the  excitation  pulse,  and  the  transmit  and  receive  aperture 

element  impulse  responses.  If  PTRz  is  the  goal  psf  at  range  z ,  the  SSE  between  the  goal  and  actual  pulse-echo  psfs  can 
be  expressed  using  equation  12,  and  minimized  to  obtain  the  optimum  receive  weights  as  shown  below. 

SSE  =  [pTR:  -  PTRz )  [pTRz  -  P TRz )  (13) 

P  ~  {Azz  Azz )  Azz  PTRz  -  Azz  PTRz  (14) 

where  aJ*  is  the  pseudoinverse  of  Azz . 

2.5  Modified  broadband  formulation  for  reduced  computational  complexity 

The  calculation  of  aperture  weights  in  the  MSSE  technique  requires  significant  computational  resources,  due  to  the 
pseudo-inverse  operation  and  the  large  matrices  involved.  However,  the  lateral  symmetry  of  the  apertures  and  the  psfs 
can  be  exploited  in  order  to  reduce  the  computational  complexity  of  the  broadband  formulation.  We  can,  if  we  choose, 
use  just  half  of  both  the  goal  and  actual  psfs  for  the  calculation  of  the  optimal  weightings.  The  symmetry  of  the  transmit 
and  receive  apertures  can  also  be  used  to  reduce  the  size  of  the  matrices  involved.  Pairs  of  elements  can  be  considered  by 
grouping  elements  that  are  on  opposite  sides  and  at  the  same  distance  from  the  center  axis.  The  computational 
complexity  is  therefore  reduced  by  a  factor  of  4.  This  concept  is  illustrated  in  figure  2. 

2.6  Application  to  enhance  Depth  of  Field  (DOF) 

The  Minimum  Sum  Squared  Error  (MSSE)  technique  that  is  described  above  for  CW  and  broadband  systems  is 
extremely  general  and  can  be  applied  in  wide-ranging  design  scenarios.  As  an  example,  we  describe  the  application  of 
the  technique  for  improved  Depth  of  Field. 

The  Depth  of  Field  (DOF)  of  an  ultrasound  imaging  system  is  generally  defined  as  the  axial  region  over  which  the 
system  is  in  focus,  or  more  rigorously,  the  axial  region  over  which  the  system  response  is  uniform  within  some 
predetermined  limit.  It  is  usually  desired  that  the  system  psf  remains  similar  to  the  psf  at  the  focus  for  as  large  an  axial 
range  as  possible. 

Currently  techniques  that  are  used  to  improve  the  DOF  include  apodization  and  dynamic  receive  focusing.  A  static 
apodization  function  is  generally  applied  to  the  transmit  aperture,  while  dynamic  apodization  is  implemented  on  the 
receive  aperture.  If  the  MSSE  technique  is  implemented  for  every  range  under  consideration  with  the  goal  psf  being  the 
psf  obtained  at  the  focus,  we  can  formally  derive  receive  apodization  weightings  that  force  the  psf  at  each  specific  range 
of  interrogation  to  be  maximally  similar  by  minimizing  the  SSE.  These  weightings  can  then  be  used  to  implement 
dynamic  apodization  and  maximize  the  DOF. 


3.  SIMULATIONS 


We  have  implemented  two  sets  of  simulations  in  order  to  illustrate  the  working  of  the  MSSE  technique.  The  first  set  was 
intended  to  illustrate  the  implementation  of  the  MSSE  technique  to  obtain  predetermined  system  psfs  in  CW  and 
broadband  systems.  The  second  set  was  designed  to  implement  the  technique  to  improve  the  system  DOF  in  CW  and 
broadband  simulations.  The  default  system  parameters  are  described  in  table  1.  Unless  otherwise  mentioned,  these 
parameters  were  used  in  all  simulations.  All  simulations  were  performed  in  Matlab.  Field  II,  an  ultrasound  simulation 
package  developed  by  Jensen  [5],  was  used  for  all  broadband  simulations. 


Parameter 

Value 

Number  of  elements 

32 

Element  pitch 

135  fl  m 

Focus 

1.3  cm 

Lateral  window  over  which  the  psf  was  calculated 

70  mm 

psf  window  sampling  interval 

20  flm 

Ultrasonic  wave  propagation  speed 

1540  m/s 

Center  frequency 

10  MHz 

Bandwidth  (in  broadband  simulations) 

75% 

Temporal  sampling  of  psf  (in  broadband  simulations) 

84  MHz 

Temporal  spacing  of  weights  for  each  element 
(in  broadband  simulations) 

36  ns 

Table  1:  List  of  default  parameters  used  in  the  simulations,  unless  otherwise  mentioned. 


3.1  One-way  CW  design  example 

For  all  CW  simulations,  we  adapted  the  propagation  function  from  the  Rayleigh-Sommerfeld  equation  described  in  [6]. 
The  elements  were  treated  as  point  sources.  The  parameters  that  were  used  are  described  in  table  1.  The  goal  psf  was 
chosen  to  be  a  Hann  window  of  width  7  mm.  The  optimum  transmit  weights  were  computed  using  equation  4,  and  the 
system  psf  was  calculated  after  application  of  these  weights.  Figures  3(a)  and  3(b)  show  the  goal  psf  and  the  resulting 
system  psf  respectively.  No  transmit  apodization,  other  than  the  calculated  weights,  was  used.  Figures  3(c)  and  3(d) 
display  the  calculated  transmit  aperture  weights,  and  the  magnitude  of  the  error  between  the  goal  and  resulting  psfs 
respectively.  The  error  is  the  difference  between  the  goal  and  the  obtained  psfs. 

3.2  Two-way  CW  design  example 

A  Hann  window  of  width  1  mm  was  chosen  to  be  the  goal  two-way  psf.  The  MSSE  technique  was  implemented  by 
calculating  and  applying  the  receive  aperture  weights  that  were  obtained  using  equation  8.  No  apodization  was  applied  to 
the  transmit  aperture.  Figures  4(a)  and  4(b)  depict  the  goal  psf  and  the  psf  generated  by  using  the  calculated  weights 
respectively.  Figures  4(c)  and  4(d)  display  the  calculated  receive  aperture  weights,  and  the  magnitude  of  the  error 
between  the  goal  and  resulting  psfs  respectively. 

3.3  One-way  broadband  design  example 

All  broadband  simulations  took  advantage  of  symmetry  to  reduce  computational  complexity,  as  previously  described. 
The  propagation  matrix  Az  was  constructed  using  dual  element  impulse  responses.  We  used  apodization  to  transmit  only 
on  selected  pairs  of  elements,  thus  obtaining  their  contributions  at  each  field  point,  at  each  time  point  under 
consideration.  These  were  then  used  to  form  the  propagation  function.  The  ideal  psf  was  generated  by  axially  weighting 
a  sinusoidal  signal  by  a  Hann  window,  and  multiplying  the  result  by  a  lateral  Hann  window. 

In  all  broadband  simulations,  we  downsampled  the  propagation  matrix  by  a  factor  of  3.  This  had  the  effect  of  reducing 
the  upper  cut-off  frequency  in  the  frequency  response  of  the  FIR  filter  constructed  using  the  weights  from  the 


Figure  3.  CW  one-way  design  example  showing  (a)  the  goal  psf  as  a  function  of  lateral  position,  (b)  the  system  psf  obtained  with  the 
MSSE  technique  as  a  function  of  lateral  position,  (a)  the  calculated  transmit  weights  as  a  function  of  element  number,  and  (b)  the  error 
between  the  goal  and  the  MSSE  technique  psfs  as  a  function  of  lateral  position. 


temporal  sampling  rate  of  the  psf  (84  MHz)  to  28  MHz.  This  rate  still  provided  adequate  sampling  since  the  input  pulse 
had  a  center  frequency  of  1 0  MHz  with  a  bandwidth  of  75  %. 

Using  equation  1 1,  we  calculated  the  optimum  transmit  weights  and  the  resulting  system  psf.  Figures  5(a)  and  5(b)  show 
the  goal  and  the  MSSE  technique  psfs  respectively,  both  as  a  function  of  lateral  position  and  time.  Figures  5(c)  and  5(d) 
depict  the  calculated  transmit  weights  as  a  function  of  the  element  number  and  time,  and  the  magnitude  of  the  error  as  a 
function  of  lateral  position  and  time  respectively.  The  error  was  calculated  by  computing  the  difference  between  the  goal 
and  obtained  psfs.  We  also  envelope  detected  and  peak  detected  the  psfs  to  generate  beam  profiles.  Figures  5(e)  and  5(f) 
show  the  goal  and  system  transmit  beam  profiles  respectively. 

3.4  Two-way  broadband  design  example 

We  used  the  same  goal  psf  as  in  the  one-way  example,  except  for  a  scaling  factor  that  accounted  for  the  reduction  in 
magnitude  due  to  two-way  propagation.  The  propagation  matrix  A„ ,  however,  was  different  from  the  one-way  case.  We 

generated  the  propagation  functions  by  using  the  entire  transmit  aperture  and  receiving  only  on  selected  pairs  of 
elements.  No  transmit  apodization  was  used.  All  other  parameters  were  consistent  with  table  1 . 
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Figure  4.  CW  two-way  design  example  showing  (a)  the  goal  psf  as  a  function  of  lateral  position,  (b)  the  system  psf  obtained  with  the 
MSSE  technique  as  a  function  of  lateral  position,  (a)  the  calculated  receive  weights  as  a  function  of  element  number,  and  (b)  the  error 
between  the  goal  and  the  MSSE  technique  psfs  as  a  function  of  lateral  position. 

We  then  downsampled  the  resulting  propagation  matrix  and  used  it  calculate  the  receive  weights  that  minimize  the  SSE. 
Figures  6(a)  and  6(b)  display  the  goal  psf  and  the  system  psf  obtained  with  the  MSSE  technique  respectively.  Figures  6(c)  and 
6(d)  depict  the  calculated  receive  weights  and  the  psf  error  magnitude  respectively.  As  in  the  one-way  case,  we  envelope 
detected  and  peak  detected  the  psfs.  The  goal  and  system  psf  profiles  are  displayed  in  figures  6(e)  and  6(f)  respectively. 

3.5  Enhanced  Depth  of  Field  (DOF) 

We  implemented  the  MSSE  technique  to  maximize  the  DOF  in  both  CW  and  broadband  simulations.  We  generated  the 
goal  psf,  which  was  the  psf  at  the  focus,  using  16  element  transmit  and  receive  apertures.  We  used  a  16  element  transmit 
and  32  element  receive  aperture  for  the  actual  system  psf.  i.e.  weights  were  calculated  for  32  receive  elements.  The 
system  was  focused  at  6.5  mm.  We  performed  CW  simulations  to  maximize  the  DOF  over  an  axial  window  from  0.1 
mm  to  50  mm  that  was  sampled  every  100  ju  m.  The  lateral  window  over  which  the  CW  psf  was  calculated  was  sampled 
at  5  ju  m,  for  more  accurate  computation  of  the  Full  Width  Half  Maximum  (FWHM)  of  the  mainlobe.  We  implemented 
broadband  simulations  over  an  axial  window  from  0.5  mm  to  32.5  mm  that  was  sampled  every  2  mm.  Current  ultrasound 
systems  attempt  to  improve  the  DOF  by  using  dynamic  apodization  and  dynamic  receive  focusing,  and  we  used  these  in 


Figure  5:  Broadband  one-way  design  example  showing  (a)  the  goal  psf,  and  (b)  the  system  psf  obtained  with  the  MSSE  technique, 
both  as  a  function  of  lateral  position  and  time.  The  calculated  transmit  aperture  weights  are  displayed  in  (c)  as  a  function  of  element 
number  and  time,  and  (d)  shows  the  error  between  the  goal  and  the  MSSE  technique  psfs  as  a  function  of  lateral  position  and  time. 
Subplots  (e)  and  (f)  display  the  envelope  and  peak  detected  goal  and  system  psfs  respectively,  both  as  a  function  of  lateral  position. 
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Figure  6:  Broadband  two-way  design  example  showing  (a)  the  goal  psf,  and  (b)  the  system  psf  obtained  with  the  MSSE  technique, 
both  as  a  function  of  lateral  position  and  time.  The  calculated  receive  aperture  weights  are  displayed  in  (c)  as  a  function  of  element 
number  and  time,  and  (d)  shows  the  error  between  the  goal  and  the  MSSE  technique  psfs  as  a  function  of  lateral  position  and  time. 
Subplots  (e)  and  (f)  display  the  envelope  and  peak  detected  goal  and  system  psfs  respectively,  both  as  a  function  of  lateral  position. 
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Figure  7:  CW  enhanced  DOF  simulation  results  showing  (a)  CW  correlation  curves,  and  (b)  CW  mainlobe  FWHM  curves  in  the  top 
panel.  The  middle  panel  shows  (c)  broadband  correlation  curves,  and  (d)  broadband  mainlobe  FWHM  curves.  The  lower  panel  shows 
(e)  the  image  obtained  using  control  psfs,  and  (0  the  image  obtained  using  the  psfs  generated  after  application  of  the  MSSE  technique. 
The  images  in  (e)  and  (f)  were  formed  from  lateral  psfs  at  multiple  ranges.  Each  column  consists  of  the  psf  at  a  single  range. 


control  simulations  to  establish  a  basis  for  comparison  with  the  results  obtained  using  the  MSSE  technique.  A  Hann 
window  was  used  as  a  transmit  apodization  function  in  all  simulations,  and  as  receive  apodization  during  the  generation 
of  the  goal  psf.  The  propagation  matrix  was  downsampled  to  limit  the  frequency  response  of  the  FIR  filter  constructed 
using  the  calculated  weights. 


5.  DISCUSSION 

Figure  3  demonstrates  the  use  of  the  MSSE  technique  in  the  most  basic  ultrasound  system  configuration  that  was 
simulated  i.e.  the  one-way  CW  system.  As  shown  in  figures  3(a)  and  3(b),  the  system  psf  closely  approximates  the  goal 
psf.  The  mean  magnitude  of  the  error  between  the  psfs  in  figure  3(d)  was  approximately  0.02%  of  the  mean  goal  psf 
amplitude. 

The  results  of  implementing  the  MSSE  technique  in  the  design  of  two-way  system  responses  is  shown  in  figure  4.  It  can 
be  seen  that  the  system  psf  that  was  obtained  after  the  application  of  the  MSSE  technique  and  the  goal  psf  are  quite 
similar.  The  resulting  mean  error  magnitude  shown  in  figure  4(b)  was  5.1%  of  the  mean  goal  psf  magnitude.  The  error 
was  much  worse  than  in  the  one-way  simulation  because  the  goal  psf  was  much  narrower  and  therefore  more  difficult  to 
generate. 

Figures  5  illustrates  the  results  obtained  when  the  MSSE  algorithm  was  implemented  in  one-way  broadband  simulations. 
Observation  of  figures  5(a),  5(b),  5(e),  and  5(f)  reveals  that  the  resulting  psf  has  a  good  qualitative  similarity  to  the  goal 
one-way  psf.  The  mean  error  magnitude  was  17%  of  the  mean  psf  magnitude. 

Results  from  the  two-way  broadband  simulation  are  shown  in  figure  6.  The  mean  error  magnitude  was  19.8%  of  the 
mean  goal  psf  magnitude. 

The  errors  in  the  broadband  simulations  are  quite  large,  but  it  is  worth  noting  from  figures  5(a)  and  6(a)  that  the  goal  psfs 
used  are  difficult  to  realize  using  spherical  waves.  They  would  have  been  easy  to  generate  using  plane  waves,  but  plane 
waves  are  an  unrealistic  model  of  the  ultrasonic  field  emitted  by  transducer  elements.  Field  II  uses  spherical  waves  to 
form  realistic  element  responses.  In  spite  of  the  very  challenging  goal  psf  used  here,  there  is  a  very  good  qualitative 
agreement  between  the  goal  and  system  psfs. 

The  effect  of  the  MSSE  technique  in  improving  the  Depth  of  Field  (DOF)  is  illustrated  in  figure  7.  The  DOF  was  defined 
in  terms  of  correlation  coefficients  as  the  axial  region  over  which  the  coefficients  remained  over  0.95.  In  the  CW  case, 
the  DOF  in  the  control  case  was  13  mm,  while  it  increased  by  285%  to  almost  the  entire  interrogated  range  of  50  mm 
when  the  MSSE  technique  was  applied.  In  the  broadband  simulations,  the  DOF  increased  from  8.3  mm  in  control 
simulations  to  26.7  mm  upon  the  application  of  the  MSSE  technique,  an  increase  of  222%.  The  improvement  in  the  DOF 
in  CW  simulations  can  be  clearly  seen  in  figure  7(a),  and  in  broadband  simulations  in  figure  7(c). 

The  DOF  was  also  defined  in  terms  of  the  Full  Width  Half  Maximum  (FWHM)  of  the  mainlobe.  Here  we  considered  the 
DOF  to  be  the  region  within  which  the  FWHM  stayed  within  25%  of  its  value  at  the  focus.  The  DOF  calculated  using 
the  FWHM  criterion  increased  from  1 1  mm  in  control  simulations,  to  almost  the  entire  range  of  interrogation  i.e.  50  mm. 
This  represents  a  355%  increase  in  the  DOF  when  the  MSSE  technique  was  applied.  The  CW  control  case  and  MSSE 
technique  case  FWHM  results  can  be  seen  in  figure  7(b).  Figure  7(d)  displays  the  FWHM  information  for  broadband 
simulations.  In  the  broadband  case,  the  DOF  evaluated  using  the  FWHM  increased  from  12  mm  to  around  27  mm  upon 
application  of  the  MSSE  technique,  an  increase  of  around  125%.  Therefore,  a  significant  improvement  in  the  DOF  was 
obtained  in  both  simulations. 

A  more  qualitative  assessment  of  the  efficacy  of  the  MSSE  technique  can  be  made  using  figures  7(e)  and  7(1).  Figure 
7(e)  shows  the  CW  psfs  obtained  in  control  simulations,  and  7(f)  displays  the  psfs  generated  in  CW  simulations  that 
implemented  the  MSSE  technique.  The  psfs  are  displayed  as  a  function  of  range  and  lateral  position.  It  can  be  seen  that 
there  is  significant  broadening  of  the  psf  mainlobe  with  range  in  the  control  simulation.  Figure  7(f)  clearly  demonstrates 
the  dramatic  improvement  in  the  DOF  obtained  using  the  MSSE  beamforming  technique. 


While  the  MSSE  design  method  worked  exceedingly  well  over  the  range  of  conditions  considered  here,  we  must 
exercise  some  caution  in  interpreting  these  results  since  we  only  observed  the  performance  of  the  technique  in  a  limited 
spatial  window.  We  cannot  predict  with  certainty  what  will  happen  outside  this  design  window.  Another  concern  is  the 
effect  of  assuming  a  wrong  propagation  model  in  the  derivation  of  the  optimum  aperture  weights.  Errors  such  as  a 
mismatch  in  the  assumed  and  actual  wave  propagation  speeds  may  have  an  adverse  effect  on  the  design  method.  We  are 
currently  investigating  these  and  other  similar  concerns. 

6.  CONCLUSIONS 

The  Minimum  Sum  Squared  Error  (MSSE)  technique  is  a  general  beamforming  method  that  can  be  used  to  design 
apertures  for  specific  applications.  It  enables  the  design  of  arbitrary  beam  profiles  by  calculating  the  appropriate 
optimum  aperture  weightings.  The  system  performance  is  optimized  because  the  calculated  weightings  minimize  the 
error  between  the  actual  and  desired  system  responses.  The  algorithm  can  be  readily  implemented  in  both  Continuous 
Wave  and  broadband  systems.  In  CW  systems,  the  receive  weights  can  be  implemented  by  way  of  apodization  and  time 
delays,  or  complex  weights.  In  broadband  systems,  implementation  is  analogous  to  applying  a  dynamic  FIR  filter  to  each 
channel. 

The  MSSE  technique  has  been  shown  to  be  effective  in  designing  ultrasound  systems  that  generate  arbitrary  desired 
system  responses.  Simulation  results  in  one-way  and  two-way  CW  and  broadband  systems  demonstrate  that  it  is  easy  to 
implement  and  can  be  applied  in  a  wide  range  of  potential  applications.  These  simulations  indicate  that  the  MSSE 
technique  compares  favorably  with  current  techniques  used  in  conventional  beamforming,  and  has  the  potential  to  be 
applied  in  a  variety  of  ultrasound  system  design  problems. 
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Abstract  -  The  design  of  transmit  and  receive 
aperture  weightings  is  a  critical  step  in  the 
development  of  ultrasound  imaging  systems.  Current 
design  methods  are  generally  iterative,  and 
consequently  time-consuming  and  inexact.  We  have 
previously  described  a  general  ultrasound 
beamformer  design  method,  the  minimum  sum 
squared  error  (MSSE)  technique,  that  addresses  these 
issues.  We  provide  a  brief  review  of  the  design 
method,  and  present  results  of  simulations  that 
investigate  the  performance  of  the  technique.  We  also 
provide  an  example  of  application  by  applying  the 
technique  to  improve  the  depth  of  field  in  CW  and 
broadband  ultrasound  systems. 

L  INTRODUCTION 

We  have  previously  introduced  the  minimum  sum 
squared  error  (MSSE)  beamforming  technique  [1], 
[2].  The  MSSE  technique  can  be  applied  in  arbitrary 
system  geometries  to  design  apertures  that  optimize 
beam  parameters.  It  utilizes  a  linear  algebra 
formulation  of  the  sum  squared  error  (SSE)  between 
the  system  point  spread  function  (psf)  and  the  desired 
or  goal  psf.  Minimization  of  the  SSE  yields  unique 
aperture  weights  that  maximize  the  system  psFs 
resemblance  to  the  desired  psf.  We  first  provide  a 
brief  review  of  the  technique  and  present  a  simple 
broadband  design  example.  We  also  present  the 
results  of  simulations  that  investigate  the  performance 
of  die  technique.  Finally,  we  present  an  example  of 
application  that  shows  die  ease  of  implementation  of 
the  technique  to  solve  common  design  problems. 

B.  REVIEW  OF  THE  MSSE  TECHNIQUE 

The  phase  and  magnitude  of  the  ultrasonic  field  at 
a  point  in  space  generated  by  an  ultrasound 
transducer  element  depend  upon  several  factors 
including  die  Euclidean  distance  between  the  point 
and  die  element,  the  orientation  of  die  element 
relative  to  the  point,  the  frequency  of  the  emitted 
wave,  and  frequency  dependent  attenuation  of  the 


medium.  The  field  can  be  expressed  as  a  function  of 
the  aperture  weighting  (possibly  complex),  and  a 
propagation  function  that  includes  the  effects  of  the 
above  factors.  Therefore,  the  one-way  transmit  lateral 
psf  at  the  range  z  can  be  represented  as, 


*1.1 

*2,1  •  •  ' 

*1.2 

*2.2  •  ’  • 

/I  M 

«  •  •  • 

where  S2  is  an  M  x  N  matrix  of  propagation 
functions  with  s,j  denoting  the  propagation  function 
that  determines  the  field  at  the  point  j  due  to  the  i  “ 
element,  T  is  an  N  x  1  vector  of  aperture 
weightings,  and  the  psf  PTt  is  an  M  x  1  vector. 

Let  4  represent  the  desired  one-way  psf  for  the 
application  of  interest.  We  can  then  characterize  the 
degree  of  similarity  between  PTz  and  the  actual 
system  psf,  Pn>  by  the  SSE  between  them. 
Minimizing  this  SSE  would  yield  a  system  psf 
optimally  similar  to  the  goal  psf.  Therefore, 
beamformer  design  is  simply  the  selection  of  transmit 
aperture  weightings  such  that  the  SSE  between  die 
desired  and  actual  system  psfs  is  minimized.  Using 
equation  1,  the  SSE  can  be  expressed  as  follows. 

«»-(4-4)#(4-4) 

=(s,7’-4)"(s2r-4)  (2) 

where  the  superscript  UH"  denotes  a  conjugate 
transpose  operation.  From  [3],  the  least  squares 
solution  for  the  optimal  transmit  aperture  weights  is, 

T  -  (s  *  sz )"’  S?  4  ■  s**  4  <3> 
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Parameter 

Value 

Number  of  elements 

32 

Element  pitch 

135  //m 

Focus 

1.3  cm 

Lateral  window  over  which  the  psf  was 
calculated 

70  mm/ 180° 

psf  window  sampling  interval 

20  nm! 0.01° 

Ultrasonic  wave  propagation  speed 

1540  m/s 

Center  frequency 

10  MHz 

-6  dB  Bandwidth  (broadband  simulations) 

75% 

Temporal  sampling  of  psf  (in  broadband 
simulations) 

120  MHz 

Temporal  spacing  of  weights  for  each 
element  (in  broadband  simulations) 

25  ns 

Table  1.  Parameters  used  in  simulations.  The  Cartesian 
coordinate  system  was  used  in  the  broadband  design 
example  and  the  depth  of  field  simulations.  The  polar 
coordinate  system  was  used  in  the  simulations  involving 
the  effects  of  the  design  window  and  the  wave  propagation 
velocity  error  simulations. 


where  the  superscripts  “-1”  and  “#”  denote  a  matrix 
inverse  and  a  pseudoinverse  operation  respectively. 

This  formulation  can  be  extended  to  CW  two- 
way  (receive)  aperture  design,  and  also  to  broadband 
one  and  two-way  beamforming  [1],  [2], 

in.  SIMULATIONS 

Unless  otherwise  mentioned,  the  parameters 
described  in  table  1  were  used  in  all  simulations 
described  in  this  paper.  In  CW  simulations,  the 
propagation  functions  were  derived  from  the 
Rayleigh-Sommerfeld  diffraction  equation  [4],  Field 
II,  an  ultrasound  simulation  package  developed  by 
Jensen  [5],  was  used  in  all  broadband  simulations. 

We  implemented  the  MSSE  technique  in  a  simple 
broadband  design  example.  The  goal  psf  was 
generated  by  axially  weighting  a  sinusoidal  signal  by 
a  Harm  window,  and  multiplying  the  result  by  a 
lateral  Hann  window.  Figures  1(a)  and  1(b)  show  the 
goal  and  the  generated  psfs  respectively.  Figures  1(c) 
and  1(d)  display  the  calculated  transmit  weights  and 
the  magnitude  of  the  error  respectively.  The  error  is 
the  difference  between  the  goal  and  system  psfs.  Note 
that  geometric  focal  delays  are  not  included  in  the 
weights,  and  must  be  applied  separately. 


The  MSSE  algorithm  optimizes  the  system  psf 
only  within  the  design  window.  Effects  occurring 
outside  this  window  are  ignored,  introducing  artifacts 
in  the  ultrasonic  field  generated  outside  the  window. 
We  performed  simulations  to  investigate  the  effect  of 
the  design  window  size  on  the  generated  field.  The 
goal  psf  was  a  6°  wide  Hann  window.  The  MSSE 
technique  was  implemented  for  design  window  sizes 
of  ±15°,  ±30°,  and  ±45°.  The  system  psf  was 
computed  over  a  ±  90°  window  using  the  calculated 
weights.  Figure  2  shows  the  goal  and  generated  psfs. 
The  first  column  depicts  the  goal  psfs  while  the 
second  column  shows  the  obtained  psfs  respectively. 
The  design  window  edges  are  shown  by  dotted  lines. 

An  incorrect  estimate  of  the  ultrasound  wave 
propagation  speed  will  degrade  the  performance  of  an 
ultrasound  system.  Since  the  MSSE  technique  uses 
dynamic  shift-variant  aperture  weights,  errors  in  the 
assumed  wave  speed  are  an  important  concern. 
Therefore,  we  implemented  simulations  in  which  the 
actual  wave  speed  (1540  m/s)  was  underestimated 
and  then  overestimated  by  25  m/s  and  50  m/s.  The 
goal  psf  was  a  6°  wide  Hann  window.  The  designed 
psf  for  each  assumed  velocity  is  shown  in  figure  3. 

In  order  to  provide  a  simple  example  of  the 
application  of  the  MSSE  technique  in  a  common 
design  problem,  we  implemented  the  two-way  MSSE 
technique  in  both  CW  and  broadband  ultrasound 
systems  to  improve  the  depth  of  field  (DOF).  The 
system  psf  should  ideally  remain  similar  to  the  psf  at 
the  focus  for  as  large  an  axial  range  as  possible  for  a 
large  DOF.  Therefore,  the  technique  was  applied  at 
every  sampled  axial  range  point  with  the  goal  psf  at 
each  range  point  being  the  psf  at  the  focus.  In  both 
CW  and  broadband  simulations,  we  generated  the 
goal  psf  using  16  element  apertures,  while  we  used  a 
16  element  transmit  and  a  32  element  receive  aperture 
for  the  actual  system  psf.  The  focus  was  placed  at  6.5 
mm.  A  Hann  window  was  used  for  transmit 
apodization,  and  for  receive  apodization  during  the 
generation  of  the  goal  psf.  We  performed  control 
simulations  that  included  dynamic  apodization  and 
dynamic  receive  focusing.  Figures  4(a)  and  4(b) 
depict  CW  and  broadband  correlation  coefficients 
calculated  by  correlating  the  psf  at  each  range  with 
the  psf  at  the  focus.  The  dotted  lines  show  the  focus. 
Figures  4(c)  and  4(d)  depict  images  constructed  using 
the  control  and  MSSE  technique  CW  psfs. 
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IV.  DISCUSSION 


V.  CONCLUSIONS 


Figure  1  shows  the  results  obtained  in  the  one¬ 
way  broadband  design  example.  It  can  be  seen  from 
figures  1(a)  and  1(b)  that  the  goal  and  system  psfs  are 
very  similar.  Note  that  this  goal  psf  is  quite 
challenging  since  it  lacks  the  normal  wavefront 
curvature.  Despite  the  flat  wavefront  of  the  goal  psf, 
we  were  still  able  to  approximate  it  well. 

The  effects  of  varying  the  design  window  size  are 
shown  in  figure  2.  It  can  be  seen  that  when  a  ±  15° 
window  was  used,  the  ultrasonic  field  outside  the 
window  had  large  grating  lobes  that  occurred  just 
outside  the  design  window  (shown  by  dotted  lines). 
However,  the  magnitude  of  the  grating  lobes 
decreased  dramatically  when  the  window  size  was 
progressively  increased  to  ±  45°.  Therefore,  it  can  be 
seen  that  the  size  of  die  design  window  significantly 
impacts  the  obtained  ultrasonic  field,  and  must  be 
carefully  chosen  to  suit  the  application.  Ideally  it 
would  always  cover  ±  90",  although  computational 
and  memory  requirements  may  limit  the  practical 
range. 

The  effect  of  errors  in  the  assumed  speed  of 
sound  are  displayed  in  figure  3.  It  can  be  seen  that  the 
designed  psf  was  not  significantly  altered  from  the 
psf  observed  with  the  correct  speed  (1540  m/s). 
Therefore,  the  MSSE  algorithm  is  stable  in  the  sense 
that  small  errors  in  the  assumed  wave  propagation 
speed  do  not  appear  to  result  in  a  significant 
degradation  of  performance. 

The  ability  of  the  MSSE  technique  to  improve  the 
depth  of  field  (DOF)  is  shown  in  figure  4.  The  DOF 
was  defined  in  terms  of  correlation  coefficient  as  the 
axial  region  over  which  the  coefficient  was  above 
0.99.  In  CW  and  broadband  simulations,  the  DOF 
increased  by  249%  and  325%  respectively  over 
control  DOF  on  application  of  the  MSSE  technique. 
A  qualitative  assessment  of  the  efficacy  of  the  MSSE 
technique  can  be  made  using  figures  4(c)  and  4(d).  It 
can  be  seen  that  applying  the  MSSE  technique 
dramatically  reduces  the  broadening  of  die  psf 
mainlobe  with  range  seen  in  the  control  simulation. 
This  improvement,  though,  comes  at  the  cost  of 
slightly  higher  sidelobe  levels.  However,  as  described 
in  [1],  we  can  use  a  weighting  function  to  selectively 
emphasize  sidelobes  in  the  MSSE  design  process. 


The  minimum  sum  squared  error  (MSSE) 
technique  has  been  shown  to  be  stable  and  useful  in 
designing  ultrasound  systems  with  arbitrary  system 
responses.  It  is  efficient,  since  there  is  no  iteration, 
and  requires  very  little  design  time.  One-way  and 
two-way  CW  and  broadband  simulations  demonstrate 
that  it  is  easy  to  implement  and  can  be  applied  to  a 
wide  range  of  applications.  The  MSSE  technique  has 
therefore  been  shown  to  have  significant  potential  to 
improve  ultrasound  beamforming  and  can  be 
implemented  in  any  ultrasound  application  in  which 
better  control  of  beam  parameters  is  desired. 
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ABSTRACT 

Adaptive  beamforming  has  been  widely  used  as  a  way  to 
improve  image  quality  in  medical  ultrasound  applications 
by  correcting  phase  and  amplitude  aberration  errors 
resulting  from  tissue  inhomogeneity.  A  less-studied 
concern  in  ultrasound  beamforming  is  the  deleterious 
contribution  of  bright  off-axis  targets.  This  paper 
describes  a  new  approach,  the  constrained  adaptive 
beamformer  (CAB),  which  builds  on  classic  array 
processing  methods  from  radar  and  sonar.  Given  a  desired 
frequency  response  for  the  mainlobe  beam,  the  CAB 
reduces  off-axis  signals  by  imposing  an  optimal  set  of 
weights  on  the  receive  aperture.  A  linearly  constrained 
adaptive  filter  dynamically  adjusts  the  aperture  weights  in 
response  to  the  incoming  data.  Initial  results  show  a  factor 
of  two  improvement  in  point  target  resolution  and  a  60% 
contrast  improvement  for  low  echogenicity  cysts.  The 
CAB  could  considerably  improve  cardiac  and  abdominal 
image  quality.  We  address  implementation  issues  and 
discuss  future  work. 


1.  INTRODUCTION 

The  ability  of  commercial  ultrasound  systems  to  image 
desired  targets  is  often  hindered  by  the  presence  of  strong 
off-axis  scattering.  Echoes  from  such  off-axis  targets 
generate  broad  clutter  which  can  overshadow  the  signal 
from  desired  targets,  greatly  reducing  image  contrast.  In 
cardiac  imaging,  the  ribs  act  as  highly  echoic  undesired 
targets,  In  the  abdomen  strong  echoes  from  the  bladder 
reduce  image  contrast.  A  method  to  reduce  side  lobe 
levels  and  suppress  clutter  would  improve  diagnostic 
imaging  in  these  situations. 

Most  adaptive  imaging  techniques  used  in 
medical  ultrasound  operate  by  correcting  phase  and 
amplitude  aberration  errors  to  improve  image  contrast  and 
resolution  [1,2].  We  introduce  a  new  approach  to  image 
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enhancement,  the  Constrained  Adaptive  Beamformer 
(CAB).  Unlike  other  adaptive  beamformers,  the  CAB 
calculates  beamformer  coefficients  to  minimize  the  impact 
of  bright  off-axis  targets,  not  to  correct  for 
inhomogeneities  in  the  propagation  path.  Adaptive 
beamforming  has  been  used  in  radar  and  sonar 
applications  to  reduce  noise  in  beam  side  lobes  [3],  but 
this  generally  is  done  using  recursive  methods  to  converge 
upon  a  single  ideal  set  of  aperture  weights  for  narrowband 
sources  in  the  aperture  far-field  [3,4].  For  diagnostic 
ultrasound,  the  ideal  aperture  weighting  changes 
constantly  because  of  the  poor  shift  invariance  in  the 
aperture  near-field.  The  CAB  therefore  calculates  new 
weights  dynamically  for  each  receive  focus. 

A  typical  beamformer  for  diagnostic  ultrasound 
receives  an  RF  line  from  each  channel  of  a  transducer 
array  and  applies  appropriate  delays  to  each  channel  to 
focus  the  signal  for  a  given  number  of  focal  ranges. 
Preset  system  apodization  is  often  used  to  weight  the  RF 
lines  coming  from  the  center  of  the  aperture  more  heavily 
than  those  from  the  edges.  Finally,  the  channels  are 
summed  and  envelope  detected  to  yield  a  B-Mode  image. 

The  CAB  begins  with  the  focal  delays  already 
applied,  but  replaces  the  system  apodization  with  an 
adaptive  set  of  aperture  weights  that  are  determined  from 
incoming  RF  lines.  The  weights  are  selected  to  reduce  the 
power  coming  from  off-axis  noise  sources,  which  can  be 
accomplished  by  modifying  a  classical  constrained  least 
mean  squares  (CLMS)  algorithm  [5]. 


2.  THEORY 

The  CLMS  problem  optimal  weights  are  found  by 

mmwT7?w  (l) 

w 

constrained  subject  to 

CTw=f  (2) 
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where  w  is  the  set  of  weights  imposed  on  the  aperture,  R  is 
the  autocorrelation  matrix  for  the  input  data,  C  is  the 
constraint  matrix,  and /is  the  vector  of  coefficients  which 
constrain  the  problem. 

The  CAB  technique  uses  the  ideal  system 
frequency  response  as  the  constraint  /  to  preserve  the 
desired  signal.  This  is  specified  using  a  finite  impulse 
response  (FIR)  filter  of  length  L.  The  aperture  weights  for 
each  range  are  calculated  from  a  window  of  input  data  L 
samples  long,  so  the  filter  length  strongly  influences 
computation  time  for  the  CAB.  Accordingly,  choosing  an 
appropriate  filter  is  a  tradeoff  between  computation  time 
and  precise  frequency  response.  All  results  presented  in 
this  paper  were  obtained  using  a  tenth-order  FIR  filter 
with  the  same  center  frequency  as  the  transducer. 

The  input  data  vector  for  each  set  of  aperture 
weight  calculations,  denoted  X>  is  a  concatenation  of  L 
values  for  the  N  input  channels  (i.e.  the  first  sample  for  all 
channels,  followed  by  the  second  sample  for  all  channels, 
etc.).  The  constraint  matrix  C  serves  as  an  index  for  the 
application  of  the  constraint  Filter  and  is  defined  as 

1  0-0 

1  0-0 
0  1-0 

C  = 

0  l  —  0 
0  0-1 
O  O  ...  1 

Each  column  of  the  M  x  L  matrix  C  has  N  rows  of  L 
nonzero  entries  corresponding  to  the  N  input  channels. 
The  dimension  M  is  thus  the  product  of  N  and  the  number 
of  filter  coefficients,  L.  This  approach  assures  that  a  plane 
wave  resulting  from  a  focused  target  is  subjected  to  the 
frequency  response  of  the  FIR  filter. 

Equations  1  and  2  can  be  solved  using  Lagrange 
multipliers  in  the  manner  described  in  [6]  to  yield  the 
following  equation  for  the  optimal  weight  vector  W: 

vv=/?-'C[CT«',cf/.  (4) 

The  MxM  matrix  R  is  approximated  by 

R  =  XTX +Sl ,  (5) 


} 


N  rows 


(3) 


where  S  scales  an  identity  matrix  to  approximately  100  dB 
below  the  mean  of  the  magnitude  of  X.  This  ensures  a 
well  conditioned  matrix  for  the  inversion  in  Equation  4. 

Finally,  the  weights  are  applied  to  the  input  data 
as  follows: 

y  =  wTX .  (6) 

Equation  6  yields  the  processed  data  for  one  range  of  the 
image.  The  CAB  replaces  the  summing  of  channel  data 
required  in  conventional  ultrasound  beamforming. 
Optimal  weights  and  the  resulting  data  are  calculated  for 
each  range  and  each  line  of  the  input  data  to  produce  the 
processed  image. 


3.  EXPERIMENTAL  METHODS 

All  experiments  were  performed  using  a  Philips  SONOS 
5500  imaging  system  operating  with  a  10  MHz  linear 
array.  Conventional  transmit  focusing  was  employed, 
though  system  apodization  was  turned  off  for  transmit.  A 
single  focal  range  (coinciding  with  the  depth  of  the 
imaged  target)  was  selected  for  receive  data.  Aperture 
growth  was  disabled  on  receive  to  maintain  constant 
aperture  size  throughout  all  ranges.  Data  was  obtained 
from  each  of  128  channels  in  succession  by  controlling 
system  apodization  and  using  custom  software  developed 
by  McKee  Poland  of  Philips  Medical  Sytems. 

For  wire  target  lateral  resolution  experiments,  a 
100  pm  steel  wire  was  imaged  in  a  water  tank  at  20°C. 
The  target  was  placed  at  a  depth  of  4  cm.  Reverberation 
in  the  tank  was  reduced  using  sheets  of  NPL  Aptflex  F28 
acoustic  absorbing  rubber  from  Precision  Acoustics,  Ltd. 

Contrast  improvement  experiments  were 
performed  using  a  Gammex  RMI  404  grayscale  tissue 
mimicking  phantom  with  graphite  scatterers  and  low 
echogenicity  cysts  approximately  4  mm  in  diameter.  To 
improve  acoustic  coupling,  a  water  standoff  was  used 
between  the  transducer  array  and  phantom  surface. 


4.  EXPERIMENTS:  POINT  RESOLUTION 

A  conventional  B-mode  image  of  a  100  pm  wire  target  is 
presented  in  Figure  la;  Figure  lb  shows  the  image  formed 
using  the  CAB.  The  control  image  is  noisy  with 
pronounced  tails  only  slightly  lower  in  intensity  than  the 
target  itself.  As  Figure  2a  demonstrates,  the  full  width 
beamwidth  at  half  maximum  (FWHM)  for  the  control  data 
is  nearly  1  mm,  with  side  lobes  impinging  on  the  main 
lobe  and  less  than  10  dB  below  the  main  lobe.  The  beam 
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Figure  2.  Lateral  beam  profiles  obtained  from  a  100  pm 
wire  at  a  depth  of  4  cm  (a)  Conventional  beamforming 
image;  (b)  CAB  image. 


Figure  1.  Images  of  100  pm  wire  targets  log  compressed 
with  a  dynamic  range  of  40  dB.  (a)  Conventionally 
beamformed  image;  (b)  CAB  image. 


profile  for  the  CAB  data  given  in  Figure  2b,  however, 
reveals  a  FWHM  beamwidth  of  only  400  pm.  The 
processed  data  also  has  easily  distinguishable  side  lobes 
which  are  suppressed  below  25  dB.  The  actual  size  and 
features  of  the  wire  target  are  much  more  clearly  depicted 
in  the  CAB  image  of  Figure  lb  than  in  the  control  image. 


5.  EXPERIMENTS:  CONTRAST  RESOLUTION 

The  magnitude  of  contrast  between  tissue  features  and 
background  speckle  has  a  profound  impact  on  the  utility 
of  medical  ultrasound  images.  By  reducing  the 
contribution  of  noise  from  directions  other  than  the  focal 
direction,  the  CAB  greatly  improves  image  contrast  for 
regions  lying  in  the  focal  direction.  Figure  3  shows 
images  of  a  4  mm  low  echogenicity  cyst  phantom.  The 
cyst  is  apparent  in  the  center  of  the  control  image  (Figure 
3a),  but  the  image  is  noisy  and  the  contrast  between 
background  and  cyst  is  poor.  The  CAB  image  of  Figure 
3b  shows  the  cyst  more  clearly  and  is  far  less  noisy.  The 
contrast  ratio  for  the  unprocessed  image,  calculated  as  the 
ratio  of  average  pixel  values  between  the  background 
speckle  and  the  cyst,  is  0.34.  The  processed  image  yields 
a  contrast  ratio  of  0.55,  a  60%  improvement  in  contrast. 
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Use  of  the  FIR  filter  alone  (without  the  CAB  algorithm) 
also  improves  cyst  contrast,  but  only  by  20%. 


(a)  (b) 

Figure  3.  Images  of  low  echogenicity  cyst  mimicking 
phantom,  log  compressed  with  a  dynamic  range  of  40  dB. 
(a)  Unprocessed  image;  (b)  CAB-processed  image. 


6.  FUTURE  WORK 

The  CAB  technique  may  have  broad  applications.  Further 
experiments  will  include  bovine  and  porcine  tissue 
imaging  and  additional  point  resolution  work.  We  are 
currently  working  to  modify  the  SONOS  5500  to  enable 
simultaneous  capture  of  all  128  beamformer  channels  over 
a  period  of  1.6  seconds.  Once  this  modification  is  in 
place,  die  CAB  will  be  applied  to  human  cardiac  and 
abdominal  imaging.  A  major  constraint  of  the  CAB  is  the 
high  computational  cost  associated  with  its  application. 
We  are  exploring  ways  to  speed  its  execution. 

The  methods  presented  in  this  paper  could  be 
effectively  paired  with  other  beamforming  and  image- 
enhancing  techniques  such  as  phase  aberration  correction 
and  angular  scatter  imaging  [7,8]. 
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7.  CONCLUSION 

Through  reduction  of  off-axis  noise,  the  Constrained 
Adaptive  Beamformer  substantially  improves  both  target 
resolution  and  image  contrast  for  wire  targets  and  low 
echogenicity  cyst  phantoms. 
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Angular  Scatter  Imaging:  Clinical  Results  and  Novel  Processing 
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Human  tissues  exhibit  variation  in 
scattering  magnitude  as  the  angle  between 
transmission  and  reception  is  changed.  These 
angular  scatter  variations  result  from  intrinsic 
acoustic  properties  and  sub-resolution 
structure.  We  have  developed  a  clinical  imaging 
system  that  uses  the  translating  apertures 
algorithm  to  obtain  statistically  reliable,  local 
angular  scatter  measurements.  The  obtained 
data  can  be  processed  to  yield  novel  images. 

A  significant  problem  with  angular 
scatter  imaging  is  limited  depth  of  field  (DOF). 
We  describe  a  new  method  to  improve  DOF  by 
applying  shift  variant  filters  to  the  data 
obtained  at  each  angle.  We  show  that  this 
approach  is  optimal  in  a  minimum  sum  squared 
error  sense.  The  filter  coefficients  used  in  the 
technique  can  be  determined  via  experiment  or 
simulation.  Unlike  prior  methods,  this 
approach  does  not  assume  a  model  for  the 
source  of  decorrelation,  rather  it  includes  all 
sources  of  decorrelation  implicitly.  We  present 
simulation  results  showing  the  improvements  in 
DOF  obtained  using  this  technique. 

We  present  experimental  angular 
scatter  data  from  phantoms  and  human  subjects. 
In  one  phantom,  designed  to  mimic 
microcalcifications  in  soft  tissue,  experimental 
data  shows  the  angular  scatter  from  500  urn 
glass  spheres  falling  off  by  50%  over  a  20 
degree  range  of  interrogation  angles.  In  the 
same  phantom  the  angular  scatter  from  50  um 
sephadex  spheres  fell  off  by  only  10%  over  the 
same  range.  In  the  human  calf  muscle, 
brightness  fell  off  by  60%  over  20  degrees, 
while  tendon  brightness  dropped  by  only  20%. 
Interestingly,  the  brightest  target  in  the 
phantom  (glass  spheres)  exhibited  the  greatest 
angular  scatter  variation,  while  the  brightest 
target  in  the  calf  (tendon)  exhibited  the  least 
angular  scatter  variation.  These  results  provide 
compelling  evidence  that  angular  scatter 
properties  are  uncorrelated  to  b-mode  image 
brightness. 
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Introduction: 

As  early  as  the  mid  1980's, 
experimental  data  indicated  that  human  tissues 
exhibit  an  intrinsic  variation  in  angular  scatter 
properties.  (Note  that  the  term  angular  scatter 
refers  to  the  variations  in  scatter  with  the  angle 
between  the  transmit  incidence  and  received 
emission,  not  anisotropic  scattering.)  While 
angular  scatter  was  extensively  explored  in  ex 
vivo  experiments,  techniques  used  at  the  time 
were  unable  to  make  meaningful  angular  scatter 
measurements  in  vivo .  We  have  recently 
described  the  use  of  the  Translating  Apertures 
Algorithm  (TAA)  for  the  acquisition  of 
spatially  localized,  statistically  robust  angular 
scatter  profiles  [1], 

We  present  initial  results  in  phantoms 
and  in  vivo,  showing  that  angular  scatter 
variations  are  significant  and  are  independent 
of  b-mode  image  brightness.  We  also  describe 
shift  variant  filters  designed  to  improve  depth 
of  field  (DOF)  in  angular  scatter  images.  We 
present  simulation  results  indicating  the 
potential  of  these  filters.  Finally,  we  discuss 
directions  for  future  work. 

Experimental  Methods: 

The  TAA  was  implemented  on  a  General 
Electric  Logiq  700MR  Ultrasound  system  by 
developing  custom  scan  software  and  employing 
a  set  of  system  software  tools  developed  at  the 
GE  Global  Research  Center.  Acquisition  angle 
was  varied  for  sequential  transmit  events  by 
simulating  system  focal  zone  updates.  The  focal 
range  was  held  steady  for  each  focal  zone  while 
the  transmit  and  receive  apodization  were 
modified  to  implement  the  TAA.  The 
experiments  presented  here  utilize  an  8 
element  active  aperture  on  both  transmit  and 
receive.  Apodization  and  dynamic  receive 
focusing  were  disabled  for  these  experiments. 
Imaging  was  performed  at  roughly  6.9  MHz 
using  a  linear  array  probe  with  205  um  element 
spacing.  An  active  aperture  of  8  elements  was 
employed  with  shifts  ranging  from  0  to  9 
elements  (each  way)  over  the  range  of 
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conditions  explored.  The  system  was  focused  at 
a  range  of  1.2  cm  and  at  ~3.0  cm  in  elevation. 

A  tissue  mimicking  phantom  was 
constructed  to  explore  the  potential  contrast  of 
angular  scatter  images.  The  phantom  consisted 
of  a  background  region  of  gelatin  with  50  um 
Sephadex  added  as  a  source  of  backscatter 
contrast.  Glass  spheres  with  500  um  diameter 
were  placed  within  the  phantom  to  mimic  the 
presence  of  microcalcifications.  The  glass 
spheres  were  suspended  by  placing  them  on  the 
interface  formed  after  one  phantom  layer  had 
hardened,  but  before  a  second  had  been  poured. 

The  myotendinous  junction  of  a  healthy 
adult  female  volunteer  was  also  imaged  to 
characterize  its  angular  scatter  properties. 


Figure  1:  B-Mode  image  and  accompanying 
angular  scatter  plot  from  glass  sphere  phantom. 
Note  that  the  brighter  glass  sphere  region  exhibits  a 
more  pronounced  reduction  in  scatter  with  angle. 


Experimental  Results: 

Results  from  the  glass  sphere  phantom 
are  shown  in  figure  1.  The  left  panel  depicts  B- 
Mode  image  and  accompanying  angular  scatter 
plot  from  glass  sphere  phantom.  Background 
region  consists  of  50  um  sephadex  particles 
within  a  gelatin  matrix.  Highly  echogenic  500 
um  glass  spheres  are  visible  to  the  left  of 
center.  B-mode  image  shows  boxed  regions  that 
were  used  to  determine  image  brightness  with 
acquisition  angle 
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Figure  2:  B-Mode  image  and  accompanying 
angular  scatter  plot  from  the  myotendinous 
junction  of  the  gastrocnemius  muscle  of  a  healthy 
female  volunteer. 

Figure  2  depicts  results  from  the 
myotendinous  junction  of  the  gastrocnemius 
muscle.  The  background  region  consists  of 
skeletal  muscle  with  the  target  region 
consisting  of  connective  tissue.  Data  was 
acquired  with  the  image  plane  perpendicular  to 
the  muscle/tendon  fibers.  These  results 
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indicate  that  angular  scatter  properties  vary 
among  soft  tissues.  Furthermore,  when 
compared  with  the  results  of  figure  1,  these 
results  clearly  indicate  that  angular  scatter 
properties  are  independent  of  backscatter. 

Our  current  experimental  system 
obtains  data  over  seven  interrogation  angles. 
While  this  data  may  be  processed  in  a  variety  of 
ways,  we  have  begun  our  investigations  by 
examining  difference-weighted  images 
computed  by  envelope  detecting  the  difference 
between  the  IQ.  data  set  obtained  at  angle  es  and 
the  IQ  data  set  obtained  at  backscatter.  An 
example  d-weighted  image  for  the  previously 
described  glass  sphere  phantom  is  shown  below 
in  figure  3.  A  control  B-Mode  image  is  provided 
for  comparison.  Both  images  are  log  compressed 
to  40  dB.  Quantitative  analysis  shows  that  the 
glass  sphere  region  appears  with  ~15  dB  greater 
contrast  in  the  d-weighted  image. 


Figure  3:  B-Mode  image  (top)  and  accompanying 
d* weighted  image  (bottom)  from  glass  sphere 
phantom.  Both  are  log  compressed  to  40  dB.  The 
d-weighted  image  exhibits  ~  15  dB  higher  contrast. 


Shift  Variant  Filtering: 

While  the  TAA  maintains  a  highly 
uniform  system  response  with  angle,  it  is 
subject  to  variations  (especially  at  ranges  far 
from  the  focus)  that  limit  depth  of  field  and 
reduce  image  quality.  As  we  show  below,  the 
impact  of  these  variations  can  be  reduced  by 
applying  a  shift  variant  filter. 

We  begin  by  considering  two  T  sample 
temporal  signals,  r j  and  r2,  received  using 

different  aperture  geometries.  These  signals 
will  be  most  useful  for  differentiating  angular 
scatter  information  if  they  are  identical  for  a 
medium  with  no  angular  scatter  variations.  We 
represent  the  signals  received  from  such 
omnidirectional  targets  as  and  r2.  We 

quantify  the  similarity  of  these  signals  using 
the  mean  squared  error  between  them. 
Representing  the  signals  r{  and  f2  as  column 

vectors,  the  expected  value  of  the  sum  squared 
error  between  them  is: 

We  consider  the  signals  to  be  the  result  of  the 
interaction  between  a  shift  variant  imaging 
system  and  a  field  of  scatterers.  This  can  be 
represented  as  the  multiplication  of  a 
propagation  matrix  and  a  scattering  vector 

r}  **  Ps 

r2-Qs 

Where  s  is  the  scattering  vector  of  X  samples 
and  P  and  Q  are  T  by  X  propagation 
matrices.  The  vector  s  contains  the  amplitude 
of  the  scattering  function  throughout  space.  The 
matrices  P  and  Q  represent  the  system 
sensitivity  to  scatterers  at  specific  locations  in 
space  as  a  function  of  time.  The  product  yields 
a  time  dependent  received  signal.  To  maximize 
the  similarity  between  the  signals  we  would 
like  to  design  a  shift-variant  filter  which  can 
be  applied  to  r2  to  minimize  the  sum  squared 
error  between  fx  and  r2.  We  represent  this 
filter  as  a  compensating  matrix  C,  of 
dimensions  T  by  7\  Thus  the  compensated  sum 
squared  error  is: 

55£c-(ft-C?J(n-Cr2)) 
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substituting  the  definitions  of  rt  and  r2  into 
this  expression  and  regrouping  terms  yields: 

SSEc-(s'(P-CQ)'(P-CQ)s) 

To  simplify  intermediate  steps  we  let 
A-{P-CQ)\P-CQ)  so  that: 

SSEC  =  [s' As) 

Changing  this  expression  from  matrix  notation 
to  summation  notation  yields: 


Since  A  is  deterministic,  the  expected  value 
operator  can  be  brought  inside  the  summation 
to  yield: 

SSEc-i$(v<)\ 

;-i  j- 1 

If  we  assume  that  the  scattering  function  s  is  a 
white  noise  process  then 

Substituting  this  result  in  and  assuming  that 
« 1  yields: 

SSEC-  $4 

M 

If  we  represent  A  as  BfB  then  the  sum  squared 
error  is: 


If  we  reshape  B  Into  a  column  vector  b  such 
that  Bkj  =  bj+j^ |)  then  the  sum  squared  error  is 

simply: 

SSEc-b’b 

Following  our  earlier  definition  of  B, 
6,+T(m)  "  Py  *  •  Thus  the  sum  squared 

A-l 

error  can  be  represented  as: 


where  cr(wW-C#  and  Pr^-Pg.  This 

equation  is  a  standard  least  squares  problem 
and  can  be  readily  solved  for  the  compensating 
weights,  expressed  as  the  vector  c. 

This  problem  can  be  further  simplified 
by  dividing  it  into  a  set  of  T  independent  least 
squares  problems.  We  begin  by  considering  the 
vectors  p  and  c  to  be  block  vectors,  each 
consisting  of  T  column  vectors  denoted  p 
through  Pj.  and  c,  through  cT  respectively. 
Using  this  formalism  the  above  expression  can 
be  rewritten  as: 


From  this  expression  it  is  clear  that  the  sum 
squared  error  between  p  and  Q'c,  does  not 
depend  upon  the  other  pn  or  the  other  c„.  Thus 
the  sum  squared  error  for  any  individual 
compensating  vector  cn  is  given  by: 

ss£,..(p.-e'c,)'(p.-e'c.) 

Determination  of  the  weights  needed  to 
minimize  this  sum  squared  error  is  a  well 
known  problem  with  the  solution  equal  to: 

K-(QQ'Yop.-(&)'  p. 

where  (Q'j  is  the  pseudoinverse  of  the 

transposed  propagation  matrix  Q' .  Using  this 
method  it  is  possible  to  determine  the  required 
filter  coefficients  to  be  used  for  each  receive 
time  of  interest.  Since  this  formalism  has 
allowed  for  shift  variant  filtering,  the  filter 
coefficients  will  need  to  be  updated  for  each 
new  output  time. 

Simulation: 

The  Field  n  program  [2]  was  used  to  test 
the  potential  utility  of  the  shift  variant 
filtering  method  described  above.  We  utilized 
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geometry  and  operating  frequency  matching 
that  used  in  our  experiments.  We  generated 
lateral-axial  system  responses  at  each  time 
range  by  resampling  the  normal  Field  space- 
time  output 


Range  (mm) 


Figure  4:  Solid  line  depicts  the  correlation  between  the 
backscatter  data  and  the  angular  interrogation  with  a  shift  of 
2  elements.  The  dashed  line  indicates  the  correlation 
resulting  from  application  of  shift  variant  filters. 
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Figure  5:  Solid  line  depicts  the  correlation  between  the 
backscatter  data  and  the  angular  interrogation  with  a 
shift  of  8  elements.  The  dashed  line  indicates  the 
correlation  resulting  from  application  of  shift  variant 
filters. 


Results  of  the  application  of  the  shift 
variant  filters  are  shown  with  control  curves  in 
figures  4  and  5.  The  application  of  shift  variant 
filters  dramatically  improves  the  depth  of  field 
of  angular  scatter  imaging  without  the  large 
computational  cost  associated  with  the  advanced 
beamforming  methods  we  have  previously 
described  [3,  4], 


Conclusion: 

The  TAA  can  be  successfully  applied  to 
the  measurement  of  tissue  and  phantom  angular 
scatter  properties.  Such  properties  provide  a 
new  source  of  image  contrast.  Shift  variant 
filters  should  improve  image  depth  of  field  and 
reduce  artifacts. 

Acknowledgements: 

We  would  like  to  thank  K.  Wayne  Rigby, 
Carl  Chalek,  Anne  Hall,  and  Kai  Thomenius  of 
GE  Medical  Systems  for  their  technical  support 
on  this  project.  We  would  also  like  to 
acknowledge  the  financial  support  of  the  Susan 
G.  Komen  Foundation  for  Breast  Cancer  Research 
and  the  US  Army  Congressionally  Directed 
Medical  Research  Program. 

REFERENCES 

[1]  W.  F.  Walker  and  M.  J.  McAllister,  "C- 
and  D-Weighted  Imaging  Using  the 
Translating  Apertures  Algorithm,"  IEEE 
Transactions  on  Ultrasonics, 
Ferroelectrics,  and  Frequency  Control, 
vol.  48,  pp.  452-461,  2001. 

[2]  J.  A.  Jensen  and  N.  B.  Svendsen, 
"Calculation  of  pressure  fields  from 
arbitrarily  shaped,  apodized,  and 
excited  ultrasound  transducers,"  IEEE 
Transactions  on  Ultrasonics, 
Ferroelectrics,  and  Frequency  Control, 
vol.  39,  pp.  262-267,  1992. 

[3]  K.  Ranganathan  and  W.  F.  Walker,  "A 
Novel  Beamformer  Design  Method  for 
Medical  Ultrasound:  Part  II 
(Simulation),"  IEEE  Transactions  on 
Ultrasonics  Ferroelectrics  &  Frequency 
Control,  accepted. 

[4]  K.  Ranganathan  and  W.  F.  Walker,  "A 
Novel  Beamformer  Design  Method  for 
Medical  Ultrasound:  Part  I  (Theory)," 
IEEE  Transactions  on  Ultrasonics 
Ferroelectrics  &  Frequency  Control, 
accepted. 


2002  IEEE  ULTRASONICS  SYMPOSIUM  — 1569 


Angular  Scatter  Imaging  in  Medical  Ultrasound 


William  F.  Walker  and  M.  Jason  McAllister 
Department  of  Biomedical  Engineering 
University  of  Virginia,  Charlottesville  VA 


Abstract 

Ultrasonic  imaging  plays  a  critical  diagnostic 
role  in  a  broad  range  of  medical  specialties ;  however, 
there  continue  to  be  areas  of  clinical  medicine  where  the 
advantages  of  ultrasound  cannot  be  brought  to  bear 
because  the  targets  of  most  interest  do  not  exhibit 
sufficient  image  contrast. 

We  have  invented  a  novel  imaging  method  that 
utilizes  modified  ultrasonic  imaging  equipment  to 
interrogate  a  previously  wiexploited  source  of  image 
contrast,  namely  ultrasonic  angular  scatter  variatiotis. 
This  technique  takes  advantage  of  the  fact  that  the 
ultrasonic  scattering  from  tissue  changes  with  the  angle 
between  the  incident  (transmitted)  ultrasonic  wave  and 
the  observed  (received)  ultrasonic  wave.  Angular  scatter 
variations  result  from  spatial  variations  in  the  intrinsic 
material  properties  of  density  and  compressibility,  as 
well  as  the  geometry  of  the  tissue  microstructure. 

We  present  experimental  results  from  tissue 
mimicking  phantoms  and  in  vivo  human  tissues 
indicating  that  angular  scatter  differentiates  targets  that 
are  indistinguishable  in  conventional  ultrasound  images. 
We  also  present  cases  where  angular  scatter  images 
improve  the  contrast  of  interesting  targets  relative  to 
conventional  images.  Early  in  vivo  results  from  the 
myotendinous  junction  of  the  human  gastrocnemius 
muscle  indicate  that  different  soft  tissues  have  different 
angular  scatter  profiles  and  that  these  profiles  can  be 
used  to  discern  between  tissues. 

1.  Background: 

1.1  Theory: 

Medical  ultrasound,  like  other  coherent  imaging 
modalities  such  as  RADAR  and  SONAR,  relies 
primarily  on  backscatter  for  image  contrast.  That  is, 
ultrasound  images  depict  the  echoes  returned  to  the 
ultrasonic  transmitter,  while  neglecting  echoes 
propagating  in  other  directions.  Although  conventional 
ultrasound  images  have  incredible  diagnostic  utility, 
there  are  numerous  clinically  relevant  targets  for  which 
ultrasonic  backscatter  fails  to  provide  adequate  image 


contrast.  In  these  cases  it  may  be  helpful  to  consider 
angular  scatter. 

Angular  scatter  refers  to  the  variations  in 
echogenicity  that  occur  as  the  angle  between  the  transmit 
and  receive  beams  is  altered.  The  standard  geometry 
used  for  discussion  of  angular  scatter  is  shown  below  in 
figure  1. 


Figure  t:  Angular  scatter  geometry.  The  angle  q  is 
known  as  the  scattering  angle.  Backscatter  occurs 
at  a  scattering  angle  of  180°. 

Angular  scatter  variations  are  intrinsic  to 
acoustic  scattering,  even  for  targets  in  the  Rayleigh 
scattering  regime.  This  occurs  because  the  acoustic  wave 
equation,  unlike  the  electromagnetic  wave  equation, 
includes  two  material  properties,  the  material  density  and 
compressibility.  The  impact  of  these  properties  can  be 
seen  clearly  in  the  acoustic  Rayleigh  scattering  equation 

HI: 

p(r,0)=A  — \k2a3($^-+  ~ -^co$<?)  (1) 

v  '  r  3  l  *-  2 P'+p  ) 

where  p(r,&)  is  the  pressure  at  a  distance  r  and  a 
scattering  angle  9  relative  to  the  target  In  this 
expression  k  is  the  wavenumber  (k  =  2xfA\  a  is  the 
scatterer  radius,  x*  and  nr,  are  the  background  and  target 
compressibilities,  p  and  pe  are  the  background  and 
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target  densities,  and  A  is  an  arbitrary  constant.  This 
expression  shows  that  a  target  consisting  of  only  a 
compressibility  difference  with  respect  to  the 
background  medium  will  exhibit  no  variation  in  angular 
scatter.  Alternatively,  a  target  varying  only  in  density 
will  exhibit  an  angular  scatter  profile  following  the 
cosine  of  the  scattering  angle. 

As  the  size  of  the  scatterer  grows  beyond  the 
Rayleigh  regime,  the  variations  in  angular  scatter 
become  much  more  significant  Known  as  the  Faran 
scattering  regime  in  acoustics,  this  realm  is  analogous  to 
Mie  scattering  in  electromagnetics.  Although  scattering 
in  this  size  range  is  extremely  difficult  to  predict 
analytically,  targets  in  this  range  include 
microcalcifications,  small  calcium  crystals  associated 
with  cancer  and  atherosclerotic  plaques.  Because  these 
targets  exhibit  densities  that  are  significantly  different 
from  soft  tissue,  we  expect  that  they  will  exhibit  strong 
variations  in  angular  scatter.  Such  variations  may  prove 
to  be  a  useful  source  of  image  contrast  as  these  clinically 
important  targets  are  generally  poorly  visualized  in 
backscatter  ultrasound  images. 

1*2  Prior  Work: 

The  potential  diagnostic  utility  of  angular 
scatter  has  been  recognized  for  some  time.  Beginning  in 
the  mid-1980s  a  number  of  researchers  performed 
experiments  with  the  goal  of  quantifying  the  intrinsic 
angular  scatter  properties  of  human  tissues.  To  support 
this  goal,  these  researchers  measured  the  average  angular 
scatter  over  a  large  area,  at  a  single  frequency  [2,  3]. 
These  systems  rotated  piston  transducers  about  a  target 
to  interrogate  different  scattering  angles.  This  approach 
was  subject  to  rapid  signal  decorrelation  with  angle,  even 
when  interrogating  omnidirectional  scatterers  [4].  This 
instability  made  the  angular  scatter  profile  measured  at 
any  given  location  highly  variable,  requiring  averaging 
over  many  samples  and  spatial  locations  to  obtain 
meaningful  profiles.  The  high  degree  of  required 
averaging,  narrow  bandwidth,  and  awkward  construction 
of  these  systems  make  them  unsuitable  for  clinical 
application. 

Recognizing  the  diagnostic  potential  of  angular 
scatter,  other  researchers  developed  clinical  systems  able 
to  image  at  a  single  scattering  angle  other  than 
backscatter  [5,  6J.  These  approaches  used  one  or  more 
phased  array  transducers  with  the  transmit  aperture 
displaced  physically  from  the  receive  aperture.  By 
applying  electronic  focusing  and  beam  steering,  these 
systems  were  able  to  interrogate  a  2-D  region  at  high 
spatial  resolution  and  with  broad  bandwidth.  Once 
configured,  the  transmit  and  receive  aperture  locations 
were  fixed,  making  it  impossible  to  acquire  angular 
scatter  data  at  more  than  two  angles  (angular  scatter  and 
backscatter)  in  any  given  location.  Acquired  angular 


scatter  images  were  displayed  next  to  corresponding  B- 
mode  images,  however  direct  comparison  was  difficult 
because  each  image  presented  a  different  speckle  pattern. 
While  these  systems  have  yielded  interesting  results, 
their  failure  to  coherently  process  data  acquired  at 
different  scattering  angles  limits  their  sensitivity  to 
angular  scatter  variations. 

We  have  developed  a  new  approach  to  angular 
scatter  imaging  that  has  the  advantages  of  both  prior 
approaches,  without  their  significant  limitations.  Our 
method  uses  a  clinical  phased  array  imaging  system  to 
acquire  angular  scatter  profiles  throughout  a  plane,  in 
real  time  [4].  By  using  broadband  excitation  and 
focusing  we  are  able  to  maintain  the  high  spatial 
resolution  typical  of  standard  B-mode  images.  The  use  of 
novel  software  and  a  phased  array  transducer  allows  us 
to  acquire  data  at  multiple  scattering  angles  without 
physically  moving  the  transducer.  The  key  to  acquiring 
stable  angular  scatter  profiles  is  the  translating  apertures 
algorithm. 

2  Methods: 

2*1  The  Translating  Apertures  Algorithm: 

The  translating  apertures  algorithm  (TAA)  is  a 
method  of  data  acquisition  that  allows  interrogation  at 
multiple  angles  without  changing  the  system’s  response 
to  an  omnidirectional  scatterer.  Such  an  approach  is 
critical  if  one  wishes  to  acquire  angular  scatter  profiles 
with  high  fidelity.  The  TAA  is  implemented  by 
maintaining  a  constant  focal  location  while  displacing 
the  transmit  and  receive  apertures  by  equal  and  opposite 
amounts.  In  its  standard  implementation  the  TAA  utilizes 
a  phased  array  so  that  all  steering,  focusing,  and  aperture 
translation  can  be  performed  electronically. 

The  performance  of  the  TAA  can  be  explained 
by  examining  figure  2.  The  first  row  depicts  data 
acquisition  using  conventional  backscatter  geometry. 
The  point  spread  function  exhibits  a  flat  phase  front, 
which  is  consistent  with  the  centered  lateral  frequency 
response  of  this  geometry,  shown  in  the  third  column. 
Traditional  angular  scatter  geometry  uses  a  shifted 
receive  aperture,  as  shown  in  the  second  row.  In  this  case 
the  phase  front  of  the  point  spread  function  is  tilted  and 
the  lateral  frequency  response  is  shifted  off  center.  This 
change  in  the  system  response  would  tend  to  overwhelm 
intrinsic  variations  in  target  angular  scatter.  Finally,  in 
the  third  row  we  see  the  TAA.  By  shifting  the  transmit 
and  receive  apertures  in  equal  and  opposite  directions  the 
TAA  maintains  a  flat  phase  response  and  a  centered 
lateral  frequency  response,  enabling  angular  scatter 
measurement  with  high  fidelity. 
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Figure  2;.  Aperture  geometries,  point  spread  functions,  and  k- 
space  representations  for  angular  scatter  data  acquisition.  The 
first  row  depicts  conventional  backscatter  acquisition,  the  second 
depicts  Iracfitional  angular  scatter  acquisition,  and  the  third 
depicts  angular  scatter  acquisition  with  the  Translating  Apertures 
Algorithm. 

2.2  The  Experimental  System: 

The  TAA  was  implemented  on  a  General 
Electric  Logiq  700MR  Ultrasound  system  by  developing 
custom  scan  software  and  employing  a  set  of  system 
software  tools  developed  at  the  GE  Global  Research 
Center.  The  acquisition  angle  was  varied  for  sequential 
transmit  events  by  simulating  system  focal  zone  updates. 
The  focal  range  was  held  steady  for  each  focal  zone 
while  the  transmit  and  receive  apodization  were  modified 
to  implement  the  TAA.  The  experiments  presented  here 
utilize  an  8  element  active  aperture  on  both  transmit  and 
receive.  Apodization  and  dynamic  receive  focusing  were 
disabled  for  these  experiments.  Imaging  was  performed 
at  roughly  6.9  MHz  using  a  linear  array  probe  with  205 
•m  element  .spacing.  An  active  aperture  of  8  elements 
was  employed  with  shifts  ranging  from  0  to  9  elements 
(each  way)  over  the  range  of  conditions  explored.  The 
system  was  focused  at  a  range  of  1.2  cm  and  at  ~3.0  cm 
in  elevation. 

23  Phantoms: 

Tissue  mimicking  phantoms  were  constructed  to 
explore  the  contrast  of  angular  scatter  images.  The  first 
phantom,  referred  to  here  as  the  glass  sphere  phantom, 
consisted  of  a  background  region  of  gelatin  with  50  urn 
Sephadex  added  as  a  source  of  backscatter  contrast 
Glass  spheres  with  500  um  diameter  were  placed  within 
the  phantom  to  mimic  the  presence  of 
microcalcifications.  The  glass  spheres  were  suspended 


by  placing  them  on  the  interface  formed  after  one 
phantom  layer  had  hardened,  but  before  a  second  had 
been  poured. 

The  second  phantom,  referred  to  here  as  the 
three  wire  phantom,  consisted  of  the  same 
Sephadex/gelatin  background  with  three  100  um 
diameter  wires  suspended  perpendicular  to  the  image 
plane.  The  wires  were  steel  wire,  nylon  monofilament, 
and  a  cotton/polyester  thread.  Each  target  was  degassed 
prior  to  phantom  construction  to  eliminate  trapped  gas. 

2*4  Human  Subjects: 

The  myotendinous  junction  of  a  healthy  adult 
female  volunteer  was  also  imaged  to  characterize  its 
angular  scatter  properties.  Images  were  obtained 
perpendicular  to  the  fiber  orientation.  This  experiment 
was  performed  under  a  human  investigations  protocol 
approved  by  both  the  university  Human  Investigations 
Committee  and  the  research  sponsor. 

2J5  D- Weighted  Images: 

A  broad  variety  of  angular  scatter  images  may 
be  formed  using  the  local  angular  scatter  profiles 
acquired  by  the  TAA.  One  of  the  simplest  image  types, 
which  we  term  difference-weighted  (d-weighted)  images, 
are  formed  by  envelope  detecting  the  complex  data 
found  by  subtracting  the  echo  data  obtained  at  some 
scattering  angle  (other  than  backscatter)  from  the  echo 
data  obtained  at  backscatter.  Prior  to  taking  the 
difference  between  data  sets  we  applied  a  range 
dependent  phase  rotation  to  compensate  for  path  length 
differences  between  the  backscatter  and  angular  scatter 
aperture  geometries.  We  expect  that  d-weighted  images 
will  highlight  targets  with  significant  variations  in 
angular  scatter,  while  suppressing  targets  with  uniform 
angular  scatter. 

3  Results: 

3.1  Angular  Scatter  Profiles: 

Figure  3  shows  a  B-mode  image  and 
accompanying  normalized  angular  scatter  profiles 
obtained  from  the  glass  sphere  phantom.  The  bright 
region  slightly  to  the  left  of  the  image  center  is  the 
location  of  a  cluster  of  glass  spheres.  Note  that  the 
angular  scatter  profile  from  this  region  falls  rapidly  as 
the’  angle  of  interrogation  decreases.  This  is  as  expected 
since  these  scatterers  have  a  significantly  different 
density  from  the  surrounding  gelatin  medium.  The 
angular  scatter  profile  from  the  corresponding  speckle 
generating  region  shows  a  much  slower  change  with 
angle.  While  the  glass  and  speckle  regions  can  be  easily 
differentiated  by  their  brightnesses  at  backscatter,  the 
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different  angular  scatter  behaviors  of  these  targets  could 
be  used  to  differentiate  them  even  if  they  had  identical 
backscatter  magnitudes. 


figure  3:  B-Mode  image  and  accompanying  angular  scatter  plots 
from  the  glass  sphere  phantom.  Note  that  the  target  region  (500 
im  glass  spheres)  exhibits  rapid  changes  In  scatter  with  angle 
whBe  the  background  (50  um  Sephadex  In  gelatin)  exhibits  little 
change  In  scattering  with  angle.  The  cosine  squared  curve 
represents  a  possible  bias  due  to  the  transducer's  limited 
angular  response. 

Figure  4  depicts  a  B-mode  image  and 
accompanying  angular  scatter  plot  from  the 
myotendendinous  junction  of  the  human  gastrocnemius 
(calf)  muscle.  We  hypothesize  that  the  different  protein 
concentrations  and  structure  of  the  muscle  and 
connective  tissues  should  cause  them  to  have  different 
angular  scatter  profiles.  The  bright  region  to  the  left  of 
center  is  a  region  of  tendon  within  the  muscle.  The 
darker  surrounding  regions  are  skeletal  muscle.  The 
angular  scatter  plots  on  the  right  clearly  indicate  different 
angular  scatter  properties  for  these  two  tissue  types. 


figure  4:  B-Mode  Image  and  accompanying  angular  scatter  plots 
from  the  in  vivo  human  gastrocnemius  muscle.  The  more 
•chogenic  tendon  region  exhibits  a  nearly  flat  angular  scatter 
profile  while  the  lees  echogente  skeletal  muscle  shows  a  rapid 
reduction  In  angular  scatter.  Note  that  this  Is  the  opposite  of  the 
behavior  seen  for  the  glass  sphere  phantom. 

3*2  D-Weighted  Images: 

Figure  4  depicts  a  b-mode  and  accompanying  d- 
weighted  image  from  the  glass  sphere  phantom.  Both 
images  were  log  compressed  to  a  40  dB  dynamic  range. 
In  the  b-mode  image  the  glass  sphere  region  is  clearly 
visible  in  the  left  side  of  the  image,  however  a 
membrane,  formed  during  phantom  construction,  is  also 
visible  over  the  lull  width  of  the  phantom.  This 


environment  may  mimic  the  presence  of  connective 
tissue  near  a  calcificed  tissue  region.  In  the  d-weighted 
image  the  contrast  of  the  glass  sphere  region  is 
significantly  improved  with  respect  to  both  the 
background  and  the  membrane  between  layers. 


figure  5:  B-Mode  {upper  panel)  and  d-welghted  (lower  panel) 
Images  of  the  glass  sphere  phantom.  The  contrast  of  the  glass 
sphere  region  is  increased  by  15  dB  in  the  d-welghted  Image 
when  compared  to  a  background  region  at  the  same  range. 

Figure  6  depicts  a  b-mode  and  accompanying  d- 
weighted  image  from  the  three  wire  phantom.  In  the  b- 
mode  image  the  contrast  of  all  three  wire  targets  is 
roughly  identical,  making  differentiation  of  these  targets 
practically  impossible.  In  the  d-weighted  image  however 
the  steel  wire  exhibits  significantly  greater  contrast 
relative  to  the  background,  the  cotton/polyester  thread 
exhibits  significantly  lower  contrast,  and  the  nylon 
monofilament  exhibits  approximately  the  same  contrast 
relative  to  the  b-mode  image.  These  results  are 
summarized  quantitatively  in  table  1.  These  results 
indicate  clearly  that  angular  scatter  differentiates  targets 
that  are  indistinguishable  in  conventional  backscatter 
images. 
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Figure  6:  B-Mode  (upper  panel)  and  d-weighted  (lower  pane!) 
images  of  the  three  wire  phantom.  The  steel,  nylon,  and 
cotton/polyester  threads  (from  left  to  right)  which  are 
indistinguishable  in  the  b-mode  image  can  be  easily 
differentiated  in  the  d-weighted  image. 
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Steel 

Nylon 

Cotton/Poly 

B-mode 

21.3  dB 

21.9  dB 

23.8  dB 

D-weighted 

27.9  dB 

22.5  dB 

17.1  dB 

Table  1:  Target  contrast  in  the  three  wire  phantom.  All  contrast 
measurements  are  relative  to  a  range  matched  background 
region.  These  results  clearly  show  that  d-weighted  images 
differentiate  targets  which  cannot  be  differentiated  In 
conventional  bacicscatter  images. 


4  Conclusions: 

The  experimental  results  presented  here  indicate 
that  angular  scatter  imaging  offers  new,  potentially 
valuable  information  in  both  human  tissues  and  in  tissue 
mimicking  phantoms.  The  techniques  employed  here 
might  be  readily  adapted  to  other  coherent  imaging 
modalities  such  as  RADAR  and  SONAR,  where  angular 
scatter  may  offer  a  new  means  of  target  classificiation. 
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1.  INTRODUCTION 

isSKESi-:="— -•-* 

lobe  levels  and  suppress  clutter  would  improve  diagnostic  imaging  m  these  situations. 

Most  adaptive  imaging  techniques  used  in  medical  ultrasound  operate  by  correcting  phase  and 

errore  to  improve  image  contest  and  resolution  [1,2].  We  introduce  a  new  approach  to  image  enhancernent  toe 

Constrained  Adaptive  Beamformer  (CAB).  Unlike  other  adaptive  beamformers,  the  path. 

coefficients  to  minimize  the  impact  of  bright  off-axis  targets,  not  to  correct  for  inhomogeneitie,  inthepropag  P 

Adaptive  beamforming  of  this  sort  has  been  used  in  radar  and  sonar  applications  to 

but  this  generally  is  done  using  recursive  methods  to  converge  upon  a  single  ideal  set  of  ape  g 

sources  h>  the  aperture  f.r-f,eld  [3,4]  For  diagnostic  dynamically 

because  of  the  poor  shift  invariance  in  the  aperture  near-field.  The  CAB  therefore  calculates  new  weig  y 

for  each  receive  focus. 

A  typical  receive  beamformer  for  diagnostic  ultrasound  (shown  in  Figure  la)  obtains  an  RF  line  " 

transducer  array  and  applies  appropriate  delays  to  each  channel  to  focus  the  signal  for  a  given  number  of  focal  ranges^ 

Preset  system  apodization  is  often  used  to  reduce  sidelobe  levels  by  weig  ting  t  e  ia®s  c  detected  to  yield  a 

the  aperture  more  heavily  than  those  from  the  edges.  Finally,  the  channels  are  summed  and  envelope  detected  yie 

B-Mode  image. 

The  CAB  (Figure  lb)  begins  with  the  focal  delays  already  applied,  but  replaces  the  system  apodization  with  an  i adaptive 
set  ofapertureweights  that  are  determined  from  incoming  RF  data.  The  weights  are  selected  to  reduce  the  power 
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Figure  1 .  Schematic  diagrams  of  (a)  conventional  delay-and-sum  receive  beamformer  and  (b)  constrained  adaptive  beamformer. 

Time  delays  Ax„  are  precalculated  for  each  sample. 

coming  from  off-axis  noise  sources,  which  can  be  accomplished  bv  modifying  a  classical  constrained  least  mean  squares 
iCLMS)  algorithm  [5]. 

2.  THEORY 

->1  Algorithm 

-  e^erai  CLMS  approach  minimizes  the  mean  output  power  in  a  signal  given  a  linear  constraint.  This  is  done 
Juuv  calculation  of  optimal  weights  found  by 
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constrained  subject  to 


CTw  =  f 

where  w  is  a  set  of  weights  imposed  on  the  input  array,  R  is  the  estimated  autocorrelation  matrix  for  the  input  data,  C  is 
the  constraint  matrix,  and/is  the  vector  of  coefficients  which  constrain  the  problem.  The  CAB  technique  uses  t  e  1  ea 
system  frequency  response  as  the  constraint/to  preserve  the  desired  signal  received  by  channels  of  an  ultrasound 
transducer  array.  This  frequency  specification  is  delineated  by  a  finite  impulse  response  (FIR)  filter  of  length  L. 

The  aperture  weights  for  each  range  r  are  calculated  for  a  window  of  input  data  L  samples  long  (r  to  r  +  L  -  1),  so  die 
filter  length  strongly  influences  computation  time  for  the  CAB.  Accordingly,  choosing  an  appropriate  filter  is  a  tradeott 
between  computation  time  and  precise  frequency  response.  All  results  presented  in  this  paper  were  obtained  using  a 
tenth-order  FIR  filter  with  the  same  center  frequency  and  bandwidth  as  the  transducer.  Shorter  filter  lengths  were  una  e 
to  effectively  approximate  the  desired  frequency  response,  and  higher  orders  filters  significantly  increased  computation 
time  without  noticeable  improvement  in  CAB  output. 

The  input  data  vector  for  each  set  of  aperture  weight  calculations,  denoted  A",  is  a  concatenation  of  L  values  for  the  N 
input  channels  (i.e.  sample  r  for  all  channels,  followed  by  sample  r  +  1  for  all  channels,  etc.).  The  constraint  matrix  C 
serves  as  an  index  for  the  application  of  the  constraint  filter.  Figure  2  shows  a  simple  case  of  the  matrix  formulation  in 
Equation  2.  Each  column  of  the  NL  x  L  matrix  C  has  N  nonzero  entries  which  align  the  aperture  weights  for  one  range 
sample  to  the  corresponding  filter  coefficient.  This  approach  ensures  that  a  plane  wave  resulting  from  a  focused  target  is 
subjected  to  the  frequency  response  of  the  FIR  filter,  while  noise  signals  arriving  incoherently  at  the  beamformer  input 
are  minimized. 

1  0  0|T  [w„' 

10  0  w2l  j  i 

0  10  w,2  1 

0  10  w22  2 

0  0  1  w13 

0  0  1  h»23 

Figure  2.  Matrix  formulation  of  CTw  -  /for  the  case  of  N  =  2  channels  and  filter  length  L  =  3.  Subscnpts  of  w  denote  channel 

number  and  filter  coefficient  number,  respectively. 

Equations  1  and  2  can  be  solved  using  Lagrange  multipliers  in  the  manner  described  in  [6]  to  yield  the  following 
equation  for  the  optimal  weight  vector  W : 


w  =  /r]c[cT/r]c]~7- 


The  NL  x  NL  matrix  R  is  approximated  by 


R  =  XTX  +  61 ,  (4 

where  8  scales  an  identity  matrix  to  a  desired  level.  This  regularization  ensures  a  well  conditioned  matrix  for  the 
inversion  in  Equation  3.  Large  values  of  8 can  reduce  output  side  lobe  levels  slightly  [7],  but  also  reduce  the  dependence 
of  R  on  the  input  samples  and  therefore  decrease  the  effectiveness  of  the  CAB.  Though  the  use  of  8  values  as  high  as  20 
dB  below  the  root  mean  squared  of  A”  does  not  noticeably  affect  the  beamformer  output,  the  results  presented  in  this 
paper  were  obtained  using  the  smallest  8  which  allowed  a  well  conditioned  matrix,  approximately  100  dB  below  the 
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mean  of  X.  For  Equation  4  to  be  a  valid  representation  of  the  autocorrelation  of  input  data,  it  is  assumed  that  the 
signals  can  be  modeled  effectively  as  zero-mean  random  processes.  " 


input 


Once  the  optimal  weights  have  been  calculated  they  are  applied  to  the  input  data,  and  the  resulting  weighted  data  is  used 
to  form  the  output  image.  The  method  used  to  combine  the  weighted  data  differentiates  the  two  beamformers  discussed 
in  this  paper:  the  Single  Iteration  (SI-CAB)  and  Multiple  Iteration  (MI-CAB)  Constrained  Adaptive  Beamformers. 

2.2  Single  Iteration  CAB 

The  SI-CAB  (presented  in  previous  work  simply  as  the  CAB  [8])  is  named  as  such  because  each  pixel  of  the  output 
image  is  determined  by  one  iteration  of  the  algorithm.  In  this  approach,  the  NL  input  samples  are  weighted  and 
simultaneously  summed  across  both  channel  and  range  through  matrix  multiplication.  The  single  sample  output  for  each 
iteration  of  the  beamformer  is  then 


y  =  w7X.  (5) 

Optimal  weights  and  the  resulting  output  are  calculated  for  each  range  and  each  line  of  the  input  data,  and  the  output 
samples  are  assembled  without  interaction  to  form  an  image.  In  this  manner,  the  SI-CAB  replaces  the  summing  of 
channel  data  required  in  conventional  ultrasound  beamforming. 

2.3  Multiple  Iteration  CAB 

If  the  weighted  samples  are  summed  solely  across  channel,  one  iteration  of  the  algorithm  results  in  a  vector  of  L  output 
samples,  each  given  by 

v 

yi  =  '2/wnl'Xnt>  (6) 

n<=] 

where  y,  is  the  /th  sample  of  the  output  vectory  for  one  iteration,  n  is  the  channel  number,  and  ‘  •  ’  indicates  scalar 
multiplication.  The  output  segments  can  then  be  retained  and  combined  by  overlap-add  after  all  iterations  of  the 
beamformer  are  complete.  Each  line  of  the  output  image  is  then  constructed  as  shown  pictorially  in  Figure  3. 
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Figure  3.  Representation  of  MI-CAB  output  for  one  full  image  line  by  overlap-add.  Each  row  represents  one  iteration  of  the 
algorithm.  M  is  the  total  number  of  ranges  in  the  input  RF  data. 

Because  every  output  pixel  depends  on  L  iterations  of  the  algorithm,  this  new  approach  is  called  the  Multiple  Iteration 
CAB,  denoted  MI-CAB  in  the  remainder  of  the  paper. 
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2.4  Frequency  Space  Considerations 

The  constraint  specified  in  Equation  2  is  designed  to  yield  a  desired  frequency  response  in  the  focal  direction. 
Exploration  of  the  Fourier  domain  characteristics  of  the  CAB  imparts  intuition  about  its  performance.  By  forcing  the 
sum  of  calculated  aperture  weights  for  each  range  to  equal  the  FIR  filter  coefficients,  the  constraint  forces  the  central 
axis  frequency  content  of  the  weights  to  match  the  Fourier  transform  of  the  FIR  filter.  Consider  the  weights  for  a  single 
iteration  of  the  CAB  before  summing,  given  by 


wn 

wl2  - 

W\L 

w2l 

w22  ••• 

™2L 

wN2  ••• 

™NL 

where  W  denotes  the  ideal  aperture  weights  reshaped  into  an  N  x  L  matrix,  and  subscripts  indicate  channel  number  and 
range/filter  coefficient  number,  respectively. 


The  center  row  of  the  2D  Fourier  transform  of  W  corresponds  to  a  signal  with  constant  lateral  frequency,  i.e.  a  plane 
wave.  Comparison  of  this  center  row  to  the  Fourier  transform  of  the  FIR  filter / shows  them  to  be  equal.  As  a  result, 
any  portion  of  the  output  signal  resulting  from  a  plane  wave  input  will  have  the  same  frequency  content  as  the  specified 
FIR  filter.  Noise  signals  arriving  from  off-axis  are  not  retained  by  the  CAB,  and  the  overall  frequency  content  of  the 
final  image  will  be  nearly  as  desired.  The  extent  to  which  this  is  true  depends  on  the  CAB  implementation  used.  Since 
the  pixels  of  an  SI-CAB  output  image  are  calculated  independently,  the  spectrum  of  the  image  is  not  guaranteed  to  be 
bandlimited.  The  MI-CAB,  however,  does  consistently  produce  a  bandlimited  (in  range)  image  due  to  the  spatial- 
averaging  characteristic  of  each  output  Y. 


3.  POINT  RESOLUTION 


3.1  Experimental  methods 

All  experiments  were  performed  using  a  Philips  SONOS  5500  imaging  system  operating  with  a  50%  fractional 
bandwidth  8.7  MHz  linear  array.  128  channels  were  used  for  transmit  with  no  system  apodization  on  either  transmit  or 
receive.  A  single  transmit  focal  range  (coinciding  with  the  depth  of  the  imaged  target)  was  selected.  Conventional 
dynamic  receive  focusing  was  employed,  though  aperture  growth  was  disabled  to  maintain  constant  aperture  size 
throughout  all  ranges.  Data  was  obtained  at  a  sampling  rate  of  40  MHz  from  each  of  128  channels  in  succession  by 
controlling  system  apodization  and  using  custom  software  developed  by  McKee  Poland  of  Philips  Medical  Sytems. 

For  point  target  lateral  resolution  experiments,  a  20  pm  steel  wire  was  imaged  in  a  water  tank  at  20°C.  The  target  was 
placed  at  a  depth  of  4  cm.  Sheets  of  NPL  Aptflex  F28  acoustic  absorbing  rubber  from  Precision  Acoustics,  Ltd.  were 
used  to  limit  reflection  from  the  sides  of  the  tank,  but  no  other  efforts  were  made  to  reduce  reverberation. 

3.2  Results 

A  conventional  B-mode  image  (obtained  by  simple  summing  across  all  128  channels)  of  a  20pm  wire  is  presented  in 
Figure  4a,  log  compressed  with  a  dynamic  range  of  40  dB.  The  background  is  noisy  due  to  reverberation,  and  the  point 
target  is  accompanied  by  pronounced  tails  only  slightly  lower  in  intensity  than  the  target  itself.  Figure  4b  shows  the 
image  formed  for  the  same  data  using  the  SI-CAB  implementation.  The  MI-CAB  is  less  effective  for  wire  target 
experiments  and  results  using  it  are  not  presented  here. 

As  Figure  4c  illustrates,  the  full  width  beamwidth  at  half  maximum  (FWHM)  for  the  control  data  is  nearly  2  mm,  with 
sidelobes  less  than  5  dB  below  the  main  lobe  and  hardly  distinguishable  from  the  main  lobe.  The  beam  profile  for  the  SI- 
CAB  data  given  in  Figure  4d,  however,  reveals  a  FWHM  beamwidth  of  less  than  500  pm.  The  processed  data  also  has 
more  easily  distinguishable  side  lobes  which  are  suppressed  to  about  9  dB  below  the  main  lobe.  The  actual  size  and 
features  of  the  wire  target  are  more  clearly  depicted  in  the  SI-CAB  image  of  Figure  4b  than  in  the  control  image.  It  is 
also  interesting  to  note  that  the  greatest  reduction  in  clutter  intensity  occurs  in  the  sidelobe  tails  of  the  point  target,  a 
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Figure  4.  Results  for  20  pm  point  target  imaged  at  a  depth  of  4  cm  using  128  channels,  (a)  Conventionally  beamformed  image; 
(b)  SI-CAB  image;  (c)  Lateral  beam  profile  for  conventional  beamforming;  (d)  Beam  profile  for  SI-CAB. 


consistent  feature  of  images  produced  with  the  SI-CAB.  The  wire  was  slightly  angled  from  perpendicular  to  the 
transducer  during  imaging,  which  likely  explains  the  asymmetry  of  the  beam  plots. 

3.3  Effect  of  fewer  input  channels 

While  most  commercial  ultrasound  systems  can  receive  echoes  on  a  maximum  of  128  channels,  certain  transducers  or 
modes  of  operation  do  not  allow  data  to  be  acquired  on  all  channels.  In  other  cases,  a  smaller  aperture  size  is  desired 
(e.g.  some  elements  obstructed  by  discontinuous  acoustic  windows  [9])  or  necessary  (e.g.  imaging  between  the  ribs  in 
cardiac  applications).  Differences  in  required  image  acquisition  time  and  hard  drive  storage  space  also  make  the  use  of 
fewer  channels  an  attractive  option  in  certain  experimental  circumstances. 

The  behavior  of  the  SI-CAB  was  investigated  for  input  channel  numbers  N  =  64  and  N  =  32.  The  same  experimental 
data  was  used  as  in  the  above  images,  with  RF  samples  from  unused  channels  discarded.  Figures  5a  and  5b  show  point 
target  images  for  standard  beamforming  and  the  SI-CAB,  respectively,  obtained  from  only  the  center  64  channels  of  the 
array.  Figures  5c  and  5d  show  the  lateral  beam  profiles  for  each  image.  Using  the  center  32  channels  for  beamforming 
yields  Figures  6a  -  6d. 
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Figure  5.  Results  for  20  pm  point  target  imaged  at  a  depth  of  4  cm  using  64  channels,  (a)  Conventionally  beamformed  image; 
(b)  SI-CAB  image;  (c)  Lateral  beam  profile  for  conventional  beamforming;  (d)  Beam  profile  for  SI-CAB. 


Comparison  of  the  images  in  Figures  4-6  clearly  shows  that  while  using  fewer  channels  degrades  the  control  images,  it 
greatly  improves  the  SI-CAB  processed  images.  Table  1  summarizes  the  lateral  beamwidths  for  standard  beamformer 
and  SI-CAB  processing  of  128,  64,  and  32  channel  data. 


Full  Width  at  Half  Maximum  [pm] 

128  channels  64  channels  32  channels 

Conventional  Beamformer 

Single  Iteration  CAB 

1570  2160  3360 

480  230  220 

Table  1 .  Lateral  maximum  amplitude  beamwidths  (in  microns)  for  processed  point  target  data. 
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Figure  6.  Results  for  20  pm  point  target  imaged  at  a  depth  of  4  cm  using  32  channels,  (a)  Conventionally  beamformed  image; 
(b)  SI-CAB  image;  (c)  Lateral  beam  profile  for  conventional  beamforming;  (d)  Beam  profile  for  SI-CAB. 

4.  CONTRAST  RESOLUTION 


4.1  Experimental  methods 

Contrast  improvement  experiments  were  performed  using  a  Gammex  RMI  404  grayscale  tissue  mimicking  phantom 
with  graphite  scatterers  and  low  echogenicity  cysts  approximately  4  mm  in  diameter.  To  improve  acoustic  coupling,  a 
small  amount  of  water  was  used  between  the  transducer  array  and  phantom  surface. 

4.2  Experiments 

The  magnitude  of  contrast  between  tissue  features  and  background  speckle  has  a  profound  impact  on  the  utility  of 
medical  ultrasound  images.  By  reducing  the  contribution  of  noise  from  directions  other  than  the  focal  direction,  the  CAB 
should  improve  image  contrast  for  regions  lying  in  the  focal  direction.  Figure  7  shows  images  of  a  4  mm  low 
echogenicity  cyst  phantom.  The  cyst  is  apparent  in  the  center  of  the  control  image  (Figure  7a),  but  the  image  is  cluttered 
and  the  contrast  between  background  and  cyst  is  poor.  The  SI-CAB  image  of  Figure  7b  shows  the  cyst  more  clearly  and 
is  far  less  noisy.  Figure  7c  displays  the  image  obtained  with  the  MI-CAB  approach,  which  exhibits  an  even  greater 
improvement  in  contrast.  Though  the  number  of  channels  used  for  processing  does  affect  the  beamformer  output,  the 
difference  is  not  as  obvious  as  for  the  point  resolution  experiments;  the  results  are  thus  omitted  from  this  paper. 
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Figure  7.  Images  of  low  echogenicity  cyst  mimicking  phantom,  log  compressed  with  a  dynamic  range  of  40  dB.  (a)  Unprocessed 

image;  (b)  SI-CAB  processed  image;  (c)  MI-CAB  processed  image. 


The  contrast  ratio  for  the  unprocessed  image,  calculated  as  the  ratio  of  average  pixel  values  between i  Ae  backgr ound 

speckle  and  the  cyst,  is  1.82.  The  SI-CAB  processed  image  yields  a  contrast  ratio  inspection 

contrast.  Processing  the  data  by  MI-CAB  produces  a  73%  improvement  m  with  a  contrast  ratio  of  3 U  ^  ‘^Sp 

may  suggest  that  the  images  in  Figures  7b  and  7c  could  be  obtained  strictly  by  applying  the  FB 
CAB  algorithm),  but  this  is  not  the  case.  Use  of  the  filter  alone  does  improve  cyst  contrast,  but  on,y  ^  ^O/oXlose 
examination  of  the  images  also  reveals  that  the  CAB  has  a  selective  effect  on  the  features  in  the  control  image, 
eliminating  some  brighf  speckles  while  better  defining  others.  We  suspect  that  speckles 

are  actually  single  bright  scattered.  Further  experiments  using  spherical  scatterers  such  as  Sephadex  beads  will  be 
performed. 


5.  DISCUSSION 


The  inmkmerttatiTns  of  the  beamformer  considered  in  this  paper  each  have  strengths  and  weaknesses.  Asmentioned 
earlier  ihe  MI-CAB  offers  a  more  stable  bandwidth  for  output  images  due  to  the  smoothing  effect  of  multiple  algorithm 
iterations.  In  fact,  this  approach  is  generally  more  robust.  We  believe  that  this  explains  its  superior  performance  m 
contrast  experiments,  where  features  of  interest  are  considerably  larger  than  the  resolution  volume  of  the  transducer.  For 
point  resolution,  however,  desired  targets  are  small  and  produce  echoes  of  much  greater  amplitude  than  then- 
surroundings.  Each  SI-CAB  output  pixel’s  relative  independence  allows  it  greater  dynamic  range  in  improvi  g 
resolution  of  such  targets. 

The  Constrained  Adaptive  Beamformer  shows  the  potential  to  dramatically  improve  image  resolution  and  contrast  in 
many  important  applications.  Unfortunately,  even  for  the  types  of  ultrasound  data  presented  in  this  paper,  the  CAB  is  not 
entirely  reliable.  The  SI-CAB  is  extremely  effective  for  improving  resolution,  but  is  also  sensitive  to  airay  imperfecttons. 
The  MI-CAB  responds  well  to  individual  element  gain  aberrations,  but  is  not  as  powerful  at  reducing  side  lobe  leve 
under  ideal  conditions  as  the  SI-CAB. 
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sample  number 


Figure  8.  Simulation  results  for  processing  of  a  zero-mean  Gaussian  plane  wave  echo.  The  left  panel  shows  time  domain  signals  for 
(a)  one  channel  of  the  input  signal;  (c)  SI-CAB  output  signal;  (e)  MI-CAB  output  signal.  The  x-axis  is  sample  number  and  the  y-axis 
indicates  amplitude.  The  right  panel  shows  the  frequency  content  of  the  signals,  where  (b),  (d),  and  (f),  correspond  to  the  Fourier 
transform  of  (a),  (c),  and  (e),  respectively.  The  x-axis  is  frequency  in  MHz,  while  the  y-axis  gives  amplitude  in  dB. 

To  illustrate  this,  we  constructed  a  zero-mean  Gaussian  plane  wave  echo  with  the  same  center  frequency  as  the  FIR  filter 
(and  therefore  the  transducer  used  to  obtain  experimental  data)  to  simulate  the  echo  from  a  focused  point  target.  The 
input  signal  was  specified  as  the  broadband  constructed  pulse  (shown  in  Figure  8a)  duplicated  across  seven  input 
channels.  Figure  8c  displays  the  time  domain  output  signal  from  beamforming  with  the  SI-CAB,  and  Figure  8e  displays 
the  MI-CAB  output.  The  right  panel  of  the  figure  shows  the  same  signals  in  the  Fourier  domain. 

The  MI-CAB  output  signal  in  Figure  8f  retains  almost  the  same  frequency  content  as  the  input  signal.  The  time  domain 
amplitude  for  Figure  8e  is  slightly  lower  than  the  input  signal,  but  this  can  be  corrected  through  simple  scaling  and  is  not 
believed  to  indicate  poor  performance.  The  SI-CAB  output  signal  shown  in  Figure  8c  is  closer  to  maintaining  the  input 
amplitude,  and  the  high  frequency  components  of  the  input  signal  (above  1 0  MHz)  are  reduced  significantly. 

For  the  next  simulation,  one  channel  of  the  input  array  was  scaled  to  10%  higher  than  the  others,  approximating  an  array 
gain  imperfection.  Figure  9a  shows  the  input  signal  for  one  of  the  unaltered  channels,  and  the  right  panel  again  presents 
Fourier  domain  signals.  In  this  trial,  the  CAB  performance  is  substantially  degraded.  The  SI-CAB  time  domain  output 
given  in  Figure  9c  is  a  strongly  modulated  version  of  the  input  signal.  The  signal  amplitude  is  orders  of  magnitude  lower 
than  that  of  the  input  signal  (note  the  different  scale  of  Figure  9c),  which  indicates  a  failure  of  the  algorithm  to  retain  the 
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Figure  9.  Simulation  results  for  processing  of  a  plane  wave  echo  with  one  aberrated  channel.  The  left  panel  shows  time  domain 
signals  for  (a)  one  unaltered  channel  of  the  input  signal;  (c)  SI-CAB  output  signal;  (e)  MI-CAB  output  signal.  The  x-axis  is  sample 
number  and  the  y-axis  indicates  amplitude.  The  right  panel  shows  the  frequency  content  of  the  signals,  where  (b),  (d),  and  (f), 
correspond  to  the  Fourier  transform  of  (a),  (c),  and  (e),  respectively.  The  x-axis  is  frequency  in  MHz,  while  the  y-axis  gives 

amplitude  in  dB. 


input  plane  wave  effectively.  The  Fourier  domain  SI-CAB  output  of  Figure  9d  shows  modulated  frequency  content  as 
well  as  amplified  high  frequency  noise.  This  result  is  certainly  problematic,  as  a  10%  difference  in  amplitude  on  a  single 
channel  of  experimental  data  could  be  expected  from  electronic  noise.  To  explain  this  undesirable  behavior,  the  method 
in  which  optimal  weights  are  calculated  must  be  revisited.  When  the  SI-CAB  receives  a  signal  with  large  amplitude  on 
one  channel,  it  minimizes  the  output  power  while  meeting  the  constraint  by  imposing  a  strong  negative  weight  on  that 
channel.  This  allows  the  weights  across  all  channels  to  sum  to  the  required  filter  coefficient,  but  does  not  allow  the 
output  sample  to  adequately  reflect  the  presence  of  a  plane  wave  signal. 

The  MI-CAB  is  more  robust  to  this  problem  through  the  smoothing  capability  of  its  multiple  iterations.  Figure  9e  shows 
that  the  time  domain  output  of  the  beamformer  retains  the  same  general  shape  as  the  input  signal.  High  frequency  noise 
is  added  (see  Figure  9f),  but  the  Fourier  domain  signal  remains  smooth  and  unmodulated  in  the  frequency  band  of 
interest.  This  result  is  encouraging,  though  future  work  with  the  CAB  will  be  directed  at  further  reducing  its  sensitivity 
to  array  imperfections  and  single  channel  noise.  Initial  results  in  which  a  constraint  is  placed  on  the  maximum 
magnitude  of  any  given  weight  have  been  promising. 
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6.  CONCLUSION 


* 


Through  reduction  of  off-axis  noise,  the  Constrained  Adaptive  Beamformer  substantially  improves  both  target  resolution 
and  image  contrast  for  wire  targets  and  low  echogenicity  cyst  phantoms.  The  Single  Iteration  CAB  performs  best  for 
point  resolution  and  has  excellent  dynamic  capability.  The  SI-CAB  is  especially  effective  for  RF  data  sets  using  fewer 
than  128  channels.  The  Multiple  Iteration  CAB  is  more  robust  than  the  SI-CAB  and  has  proved  to  be  most  appropriate 
for  enhancing  contrast  resolution. 
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